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ABSTRACT 

Shape control of spacebome antennas can provide the ability to correct for effects 
such as thermal distortion and manufacturing errors as well as control the shape of an 
antenna’s radiated beam. This thesis examines the performance of piezoceramic actuators 
in producing static deformation of a cantilever beam and analyzes the optimal placement 
of actuators to best approximate a desired deformation profile. Predictions of actuator 
effectiveness at producing curvature are developed using an Euler-Bernoulli model. An 
algorithm to determine the optimal locations and input voltages for a fixed set of actuators 
to achieve a desired deformation profile of a cantilever beam using embedded Nelder and 
Mead simplex search routines is presented and evaluated for two shape functions and 
various combinations of actuators. Experimental measurements show that the Euler- 
Bemoulli model provides a reasonable prediction of actuator performance at low input 
voltage but does not account for nonlinear behavior of the piezoceramic and the effects of 
hysteresis and transverse stresses. Further experiments demonstrate the ability of four 
piezoceramic actuators to produce an approximation of a parabolic deformation profile of 
a cantilever beam and illustrate the importance of considering these effects in determining 
the required actuator input voltages. 
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I. 



INTRODUCTION 



A BACKGROUND 

Control of the shape of long flexible structures is of interest for a variety of space- 
related applications. Many spacecraft designs include slender flexible appendages as 
components of on board antenna systems or sensors. Maintaining the desired shape of 
such an appendage can be a central requirement for obtaining optimal performance from 
its associated system. Effects such as thermal distortion, manufacturing errors, spacecraft 
acceleration and vibration can all cause deviation of the shape of a component from its 
desired profile. Active means to counter these and other effects are therefore of 
considerable interest. 

Piezoceramic actuators offer an attractive means to produce deformation of 
flexible structures. These actuators have no moving parts and require only a supplied 
electrical voltage to produce an electric field across the actuator material in order to 
fimction. Piezoceramics can be bonded to or, when properly insulated, embedded inside 
the structure they are intended to control. Applying an electric field to a piezoceramic 
actuator causes deformation of the piezoceramic material which, in turn, applies stress to 
the structure it is attached to. 

While substantial research effort has been devoted to the employment and 
placement of piezoceramic actuators for the active suppression of structural vibrations, 
considerably less attention has been focused on the use and optimal placement of 
piezoceramics for quasi-static shape control of flexible structures. Although 
piezoceramics are limited to the production of relatively small structural deformations, 
they may be more than adequate for certain applications, such as countering thermal 
effects or manufacturing errors or producing a small desired deformation profile in a 
flexible structure. 
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B. SCOPE OF THESIS 



This thesis examines the performance of piezoceramic actuators in producing static 
deformation of a thin cantilever beam and analyzes the optimal placement of actuators to 
best approximate a desired deformation profile. Models of common configurations of 
piezoceramic actuators bonded to the surface of the beam are developed from the Euler- 
Bemoulli model of Crawley and Anderson [Ref 1] to compare the effectiveness of the 
different configurations and examine the effects of actuator thickness, width and material 
properties on the ability of the actuator to produce curvature of the composite 
actuator/beam structure. Predictions of actuator performance from the Euler-Bemoulli 
model are compared to experimental measurements of the tip displacement and surface 
strain of a cantilever beam for several actuator configurations. 

An algorithm to determine the optimal locations and input voltages for a fixed set 
of actuators which minimize the error between the actual and desired shape of a cantilever 
beam under piezoceramic actuator induced deformation is presented. Optimization of the 
actuator locations and voltages is accomplished using the Nelder and Mead simplex search 
algorithm for the optimization of nonlinear functions. The MATLAB software package is 
used to implement and test the algorithm for two shape functions and various 
combinations of actuators. Capabilities and limitations of the algorithm are discussed. 
The predicted performance of four actuators placed in the optimum locations determined 
by the algorithm to approximate a parabolic beam deformation profile is compared to the 
results from a second set of experiments. 
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n. PIEZOCERAMIC MATERIAL PROPERTIES 



Piezoceramics are part of a class of materials which exhibits piezoelectric behavior. 
Such materials have the property that material strain is produced in response to an electric 
field applied to the material in line with the its "poling" direction. Similarly, strain due to 
mechanical stresses on the material will produce an electric field in the material poling 
direction due to the accumulation of electrical charge on the material surfaces normal to 
the poling axis. This accumulated charge results in a voltage between the two surfaces. 
Figure 1 shows the standard coordinate system used to describe the orientation of a 
piezoelectric material. The poling axis is normally designated as the 3-axis, with the 1- 
and 2-axes indicating the principal material directions orthogonal to the poling axis. 

Piezoelectric materials are typically orthotropic, meaning that the material 
properties are uniform along each principal material axis, but the properties along one 
principal axis may not be the same as those along another principal axis. Piezoceramics, 
however, are generally transversely isotropic, meaning that the material properties are the 
same in all directions in the plane orthogonal to the poling axis, depicted by the 1-2 plane 
in Figure 1 . 




Figure 1. Piezoceramic Material Coordinate System Conventions 
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The strain behavior of piezoelectric materials in response to applied electric fields 
and the voltage produced by mechanical stresses makes piezoelectric materials suitable for 
use as both actuators and sensors in a variety of applications. The voltage between the 
two surfaces normal to the poling axis in both actuation and sensing applications requires 
that these two surfaces be electrically insulated from one another, hence piezoelectric 
materials are generally bonded to the surface of conductive metal structures with the 
material poling axis normal to the bond. Piezoelectric materials can, however, be 
embedded in electrically insulating materials, including many composites. 

A. CONSTITUTIVE RELATION 

The axial strains in a piezoelectric material are functions of the normal stresses 
applied to the material and the electric field applied in the material poling direction. These 
strains are given by the material constitutive relation; 
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where the 3-direction is the poling direction of the piezoelectric material, the 1 and 2- 
directions are the principal material axes orthogonal to the poling direction, s, is the axial 
strain in the /-direction, and dji are the material elastic moduli, Poisson ratios and 

piezoelectric coefficients, respectively, cr, is the axial stress applied in the /-direction and <{>3 
is the electric field applied in the material poling direction. The piezoelectric coefficients 
are generally dependent on the value of the applied electric field but are typically assigned 
a constant value which provides a linear approximation of the material's piezoelectric 
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behavior over a specified operating range. The material elastic moduli are also field 
dependent and are typically specified for constant electric field (short circuit) and open 
circuit conditions along the poling axis. The moduli for constant electric field are used in 
Equation 2. 1 . Typical values of the material properties for various types of lead zirconate 
titanate (PZT), the most commonly used piezoceramic material, are listed in Table 1 . 



Material 


PZT-4 


PZT-5A 


PZT-8 


PZT-5J 


PZT-5H 


DOD Type 


I 


II 


III 


V 


VI 


Ei, E2 (x 10‘° nW) 












• Short Circuit 


8.2 


6.1 


8.7 


6.2 


6.2 


• Open Circuit 


9.9 


6.9 


9.9 


7.1 


7.1 


E3 (x 10*“ NV) 












• Short Circuit 


6.6 


5.3 


7.4 


4.4 


4.8 


• Open Circuit 


12.6 


10.6 


11.8 


n/a 


11.1 


dsi (x 10'*^ mA^) 


-122 


-171 


-97 


-220 


-274 


u (approximate) 


0.31 


0.31 


0.31 


0.31 


0.31 



Table 1. Material Properties of Selected Piezoceramic Materials [After Ref 2] 



It should be noted that the piezoelectric constitutive relation is of the same form as 
the thermal-mechanical constitutive relation, with the piezoelectric coefficients 
corresponding to the material coefficients of thermal expansion and the electric field 
applied in the poling direction corresponding to the difference in material temperature 
from its reference temperature. This fact proves useful in modeling piezoelectric material 
behavior using finite element software which includes thermal-mechanical models but not 
piezoelectric models. 
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B. ACTUATOR STRESS 



If the piezoceramic material is under a state of plane stress in the plane orthogonal 
to the poling axis, i.e. 03 = 0, then the expression for the transverse strains in Equation 2.1 
can be simplified to: 




( 2 . 2 ) 



which can be solved for the transverse normal stresses; 
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where the reduced stifihesses Qij are given by: 
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and the material Poisson ratios are related by: 
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Equations 2.1 through 2.5 apply to orthotropic piezoelectric materials and must be 
used for materials such as polyvinylidene flouride (PVDF) film which are not transversely 
isotropic in the plane orthogonal to the poling axis. Transversely isotropic materials such 
as PZT ceramics have the same properties in all directions in the 1-2 plane, i.e. 

= E2 

^12 = ^21 ( 2 . 6 ) 

^31 “ ^32 

which allows Equation 2.3 to be further simplified to: 
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(2.8) 



If the material is under a state of uniaxial stress along the 1-axis, i.e. 02 = co = 0 in 
Equation 2.1, then the expression for the axial strain along the 1-axis can be reduced to: 



El 



which can be solved for the normal stress along the 1-axis: 



(2.9) 
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(J, — £^1(^1 ^31^3) 



(2.10) 



The plot of Equation 2.10 in Figure 2 is illustrative of the performance of a 
piezoelectric material as an actuator. An unconstrained piezoelectric material will have 
normal strain ei equal to the actuation strain and will therefore produce no normal 
stress a/. In contrast, a piezoelectric material constrained to have no normal strain ei will 
produce the material's blocked normal stress: 

<^Iblocked ~ ~^1^31^3 (21 1 ) 

In general, the normal stress in a piezoelectric actuator bonded to or embedded in a 
deformable structure will be somewhere between zero and the blocked normal stress, as 
determined by the material properties and geometry of the actuator and structure. 




Figure 2. Axial Stress V5. Strain For Piezoelectric Material Under Uniaxial Stress 
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m. EULER-BERNOULLI MODEL OF BEAM DEFORMATION DUE TO 
SURFACE-BONDED PIEZOCERAMIC ACTUATORS 



Piezoceramics can be used as actuators to produce curvature of a section of a 
beam. Voltage applied to a piezoceramic actuator bonded to the surface of a beam as 
shown in Figure 3 will induce normal strain of the same sign in both the 1- and 2- 
directions. The beam will resist expansion of the piezoceramic to its full actuation strain, 
producing normal stress in the piezoceramic which will cause curvature of the composite 
actuator/beam about the x- and z-axes. The curvature about the z-axis for a configuration 
such as that shown in Figure 3 will be significantly greater than the curvature about the x- 
axis due to the difference in the moments of inertia of the respective beam cross-sections. 
The resulting stress state in the actuator and the beam can be assumed to be plane stress 
provided that the material dimensions in the 3-direction are much smaller than the 
corresponding dimensions in the 1- and 2-directions. The normal stresses in the 
piezoceramic for the plane stress condition are coupled in the 1- and 2-directions by the 
Poisson effect, as indicated by Equation 2.7. Solution of the coupled equations requires 
the use of a complex plate theory model, but a simple approximation of the curvature 
about the z-axis can be obtained using Euler-Bemoulli assumptions. 
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Euler-Bemoulli beam theory assumes that the transverse components of normal 
stress 02 and are negligible compared to the axial component oy and hence uses 
Equation 2.10 to determine the normal stress in the piezoceramic along the beam axis. 
Comparison of Equations 2.7 and 2.10 reveals that the Equation 2.10 provides a lower 
estimate of the axial stress than Equation 2.7 for a given normal strain Si and applied 
electric field (f> 3 . The Euler-Bemoulli model therefore gives a conservative approximation 
of the bending moment produced by the actuator about the beam’s r-axis for a given 
curvature of the beam about that axis. 

Crawley and Anderson [Ref 1] compared two models of beam curvature due to 
induced strain actuation by piezoceramic actuators bonded to beam structures. The first 
model, initially presented by Crawley and de Luis [Ref 2], assumed uniform axial strain in 
a pair of piezoceramic actuators bonded symmetrically to the outer surfaces of a beam 
undergoing actuator-induced bending. The second model assumed that the beam behaves 
as an Euler-Bemoulli beam, with a linear distribution of axial strain throughout the 
composite actuator-beam cross-section, as shown in Figure 4. The predicted beam 
curvature under both models was compared to that predicted by a detailed two- 
dimensional finite element model. The Euler-Bemoulli model was found to provide results 
that are within 0.1% of the finite element model over the full range of possible actuator- 
to-beam thickness ratios, while the uniform strain model results diverged significantly for 
actuator thicknesses greater than about 20% of that of the beam. The Euler-Bemoulli 
model was therefore judged to provide a much more accurate prediction of beam 
curvature produced by surface bonded piezoceramic actuators. 

Chaudry and Rogers [Ref 4] performed a similar analysis of an actuator bonded to 
one side of a beam using a thermal-mechanical finite element model incorporating two- 
dimensional linear plane stress elements to model the effect of the actuator’s piezoelectric 
behavior. The Euler-Bemoulli and finite element models were again found to be in close 
agreement, differing by at most 3.9% for actuator-to-beam thickness ratios of 1 or less. 

Bronowicki and Betros [Ref 5] refer to experimental work conducted by Griffin 
and Denoyer [Ref 6] which demonstrated that the tip displacement of a beam with a pair 
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of piezoceramic actuators mounted symmetrically on opposite sides was 8% less than 
predicted by the Euler-Bemoulli model. The difference between the predicted and 
measured displacements was within the range of measurement and manufacturing 
tolerances of the experiment. 







Figure 4. Linear Strain Distribution for Euler-Bernoulli Beam 



All of the above references support the use of Euler-Bemoulli assumptions as a 
suitable model for analysis of the effectiveness of piezoceramic actuators in producing 
curvature of the composite actuator/beam. The development of this model by Crawley 
and Anderson [Ref. 1] was used as a basis for comparison of various actuator 
configurations in this thesis. 



A. GENERAL FORMULATION 

Crawley and Anderson [Ref 1] expressed the total strain energy in an Euler- 
Bemoulli beam with attached or embedded piezoceramic actuators as; 




K 



\EA\ (£S); 

(£S). (£?),. 
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-2[Pa 




(3.1) 
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where / is the length, So is the axial strain at >^ = 0, as shown in Figure 4, and k is the 
curvature of the actuator/beam section, which can be approximated by: 



d0 

K = 

dx 



(e«i) 



(3,2) 



The total area stiffhess (EA)c, the first moment of inertia (ES)c and the second 
moment of inertia (El)c of the composite actuator/beam cross-section are respectively 
given by: 



(EAl=jE(y)«'(y)dy 

y 



(3.3) 



{ES}^ = j E(y)h'(y)yify 

y 



(3.4) 



(^)c = 1 ^(y)^(y)y^^ 

y 



(3.5) 



where E(yJ is the elastic modulus and w(yj is the width of the composite actuator beam 
material as a function of the coordinate y. The first term inside the integral in Equation 
3.1, with the 1/2 multiplication factor, represents the strain energy per unit length for a 
composite beam with the specified cross-sectional properties composed of materials which 
exhibit linear stress-strain behavior in accordance with Hooke's Law. The second term 
and 1/2 multiplication factor represent a reduction in the strain energy per unit length due 
to the portion of the material strain which is induced by the piezoelectric effect rather than 
applied stress. 

The resultant force Pa and resultant moment Ma due to the actuation strain A are 
respectively determined fi-om: 
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(3.6) 



= j E(y)My)w{y)cfy 

y 



Mj^ = ^ E{y)K(y)wiy)ydy 

y 



(3-7) 



where the actuation strain A is determined by the actuator piezoelectric charge constant 
dsi, the actuator thickness tp and the actuator input voltage V applied in the poling 
direction of the piezoceramic; 



A — 



V 



(3.8) 



If the actuator and beam properties are constant over the length of the 
actuator/beam section. Equation 3.1 can be integrated to give; 



U = Ue,\EA)^ + 2 s,k{ESI +k\EI)^ - - 2M ,k\ 



(3.9) 



Equations for the axial strain So at>^0 and the curvature k can then be obtained by 
differentiating Equation 3.9 with respect to So and fcand using Castigliano's First Theorem 






dqi 



- = Q, (/=1.2,...,«) 



(3.10) 



where qt is the /th displacement and Qi is the /th generalized external force. If no external 
forces or moments are applied to the beam, substitution of Equation 3.9 into Equation 
3.10 gives: 
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(3.11) 



= \{2s,(EA)^ + 1k(ES), - 2/\] = 0 
^ = i/[2£„(£S). +2^(£7), -2AY^]= 0 



which can be solved to obtain expressions for the axial strain at >^ = 0 and the curvature of 
the composite beam; 
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If the composite beam cross-section is symmetrical with respect to the 2 -axis, then 
the first moment of inertia of the cross-section is zero and Equations 3.12 simply to; 
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The location of the neutral axis, which has zero axial strain as shown in Figure 4, 
can be determined by solving; 

= 0 ( 3 - 14 ) 

for>/, giving; 
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(3.15) 
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B. 



MODELS OF SELECTED ACTUATOR/BEAM CONFIGURATIONS 



Piezoceramic actuators are typically bonded to beams in one of four 
configurations: a single actuator, a pair of actuators bonded symmetrically to opposite 
sides of the beam, a stacked pair of actuators bonded to one side of the beam and stacked 
pairs of actuators bonded symmetrically to opposite sides of the beam. Prepackaged 
actuators, such as the PZT QuickPacks from Active Control eXperts, Inc., which have the 
piezoceramics and wiring encased in an insulating matrix material, can also be bonded to 
one or both sides of a beam. It is of interest to compare the effectiveness of each of these 
configurations in producing curvature of the composite beam. The general formulation of 
the Euler-Bemoulli model is applied to each of these configurations in the following 
subsections. While the bond between an actuator and a beam is often negligibly thin 
when compared to the thicknesses of the beam and actuator, Baz and Poh [Ref 7] noted 
that the thickness and material properties of the bond may contribute significantly to the 
composite beam’s mechanical behavior. The dimensions and modulus of the bond material 
are therefore included in the models for completeness. For simplicity, all piezocCTamic 
actuators in the models with more than one actuator are assumed to have the same 
dimensions. 

1. Beam Deformation Due To Single Piezoceramic Actuator. 

Figure 5 shows a schematic representation of a piezoceramic actuator bonded to a 
section of a beam which extends in the ±x-direction. The total area stiffhess for the 
composite single actuator configuration is obtained by substitution of the elastic moduli Ei 
and dimensions w, and i, into Equation 3.3, where no subscript, the subscript b and the 
subscript p denote the properties of the beam, bond and piezoceramic actuator, 
respectively: 
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Figure 5. Schematic of Single Piezoceramic Actuator Bonded to Beam 



The first and second moments of inertia of the composite beam cross-section are similarly 
obtained from Equations 3.4 and 3.5: 
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Wc = j Ewy^dy + j E^w^y^dy + | E^w^dy 
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The resultant force Pa and resultant moment due to an actuator input voltage V are 
obtained from Equations 3.6, 3.7 and 3.8: 
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(3.19) 



(3.20) 



Since the composite beam cross-section is not symmetrical with respect to the z- 
axis, an actuator input voltage will produce both non-zero axial strain £b at = 0 and 
curvature k about the z-axis, as given by Equations 3.12, with the location of the neutral 
axis given by Equation 3.15. 

2. Beam Deformation Due To Symmetric Piezoceramic Actuator Pair 

Figure 6 shows a schematic representation of a symmetric piezoceramic actuator 
pair bonded to a section of a beam which extends in the ±x-direction. The total area 
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stiffness and first and second moments of inertia for a beam with a symmetric actuator pair 
can be computed as for the single actuator configuration using Equations 3.3 through 3.5; 
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Figure 6. Schematic of Symmetric Piezoceramic Actuator Pair Bonded to Beam 



The resultant force Pa and resultant moment M/i, as given by Equations 3.6, 3.7 and 3.8, 
depend on the input voltages F/ to the actuator on the +y side and F^ to the actuator on 
the -y side of the beam; 
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Since the composite actuator/beam cross-section is symmetric with respect to the 
z-axis, the axial strain at>> = 0 and curvature about the z-axis are given by Equations 3.13. 
Substitution of Equations 3.24 and 3.25 into Equations 3.13 gives: 
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(3.26) 



Thus, equal voltages applied to each actuator will produce axial strain with no curvature 
and input voltages of equal magnitude and opposite polarity will produce curvature with 
no axial strain at^ = 0. Input voltages of unequal magnitude will produce a combination 
of curvature and axial strain at_y = 0. 



3. Beam Deformation Due To Stacked Pair of Piezoceramic Actuators 

Figure 7 shows a schematic representation of a stacked piezoceramic actuator pair 
bonded to a section of a beam which extends in the ±x-direction. The total area stiffness 
and first and second moments of inertia for a beam with a stacked actuator pair are 
determined by substitution of the appropriate cross-sectional dimensions and properties 
into Equations 3.3 through 3.5; 
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Figure 7 . Schematic of Stacked Piezoceramic Actuator Pair Bonded to Beam 
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The resultant force Pa and resultant moment Ma, as given by Equations 3.6, 3.7 and 3.8, 
depend on the input voltages Vj to the inner actuator and V 2 to the outer actuator; 
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(3.31) 



Since the composite actuator/beam cross-section is not symmetric with respect to 
the z-axis, the axial strain at _y = 0 and curvature about the z-axis are given by Equations 
3.12. Actuator input voltages will, in general, produce both non-zero axial strain eo at y- 
0 and curvature x about the z-axis, with the location of the neutral axis given by Equation 
3.15. 
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4. Beam Curvature Due To A Symmetric Pair of Stacked Piezoceramic 
Actuators 

Figure 8 shows a schematic representation of a symmetric pair of stacked 
piezoceramic actuators bonded to a section of a beam which extends in the ±x-direction. 
The total area stiffiiess and first and second moments of inertia for a beam with a stacked 
actuator pair are determined by substitution of the appropriate cross-sectional dimensions 
and properties into Equations 3.3 through 3.5; 
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Figure 8. Schematic of Symmetric Pair of Stacked Piezoceramic Actuators Bonded to 

Beam 
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(3.34) 



The resultant force Pa and resultant moment Ma, as given by Equations 3.6, 3.7 and 3.8, 
depend on the input voltages Vi and V 2 to the inner actuators and V 3 and V 4 to the outer 
actuators, where the odd numbered voltages are applied to the actuators on the +y side 
and the even numbered voltages are applied to the actuators on the -y side of the beam, 
respectively: 
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(3.36) 




Since the composite actuator/beam cross-section is symmetric with respect to the 
z-axis, the axial strain at_v = 0 and curvature about the z-axis are given by Equations 3.13. 
Substitution of Equations 3.35 and 3.36 into Equations 3.13 gives: 



Hence, beam curvature with no axial strain at j = 0 can be obtained by setting Vj = -V 2 
and V3 = -V4. Similarly, pure axial strain with no beam curvature is achieved when F/ = Fj 
and Vs =¥ 4 . 

5. Beam Deformation Due To Single ACX QuickPack Actuator 

Figure 9 shows a schematic representation of a single ACX QuickPack actuator 
bonded to a section of a beam which extends in the ±x-direction. A QuickPack actuator 
consists of a stacked pair of piezoceramic wafers encased in a matrix material. 
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Figure 9. Schematic of Single ACX QuickPack Actuator Bonded to Beam 

The total area stiffness and first and second moments of inertia for a beam Avith a 
single ACX QuickPack actuator are determined by substitution of the appropriate cross- 
sectional dimensions and properties into Equations 3.3 through 3.5; 



2 2 2 
{EA)^ = j Ewdy + j E^w^dy -t | E^w^dy 

t t t 

~2 2 2 ^^* 

+ \{E^-E^)w^dy+ \(Ep-E^Wj,dy (3.38) 

= Ewt + E^wf^ + Ejvj^ + i(Ej, - E^wf^ 



27 








t 



I I 2 

{ES)^ = j + J E^^bydy + J Ejv^ydy 






r+^i+^wi+^/Fi2+2r^ 



\{^p-E„)wpydy+ j - E„ )wj,ycfy 

' " |+'4+'ml+'«2 + V 

f tA ^ ( t t. 

- + — I + E w t \ - 
2 2. 

+{Ep -E^ptp{t +2/j, +2^„, +2^p) 



= - + — I +£„w„?„| - + 



(3.39) 



t 



(£/)^ = J + j j E^wydy 



t 

r'‘ 









+ \{Ep-E„)wpy^dy 



-+^i+^nil+^fn2+^p 






12 
+£ wj 






7o 


ro 




— 


+ - 


h + — 


Vi) 


Kl) 


* 3 



m mm 



2*^ V ^V2^ 3 






“ + ^6 + ^ml J 1^2 ^ ^ ^ 



? + 2/j, 4- 2/„, + /^2 + 



^1 
3 J 



(3.40) 



The resultant force Pa and resultant moment Ma, as given by Equations 3.6, 3.7 and 3.8, 
depend on the input voltages Vi to the inner actuator and V 2 to the outer actuator: 
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Since the composite actuator/beam cross-section is not symmetric with respect to 
the z-axis, the axial strain at _v = 0 and curvature about the z-axis are given by Equations 
3.12. Actuator input voltages will, in general, produce both non-zero axial strain So 2 Ay = 
0 and curvature k about the z-axis, with the location of the neutral axis given by Equation 
3.15. 



6. Beam Deformation Due To Symmetric Pair of ACX QuickPack 
Actuators . 

Figure 9 shows a schematic representation of a symetric ACX QuickPack actuator 
pair bonded to a section of a beam which extends in the ±r-direction. The total area 
stiffiiess and first and second moments of inertia for a beam with a symmetric ACX 
QuickPack actuator pair are determined by substitution of the appropriate cross-sectional 
dimensions and properties into Equations 3.3 through 3.5: 
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Figure 10. Schematic of Symmetric ACX QuickPack Actuator Pair Bonded to Beam 
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The resultant force Pa and resultant moment Ma, as given by Equations 3.6, 3.7 and 3.8, 
depend on the input voltages Vj and V 2 to the inner actuators and V 3 and V 4 to the outer 
actuators, where the odd numbered voltages are applied to the actuators on the +y side 
and the even numbered voltages are applied to the actuators on the -y side of the beam, 
respectively: 
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Since the composite actuator/beam cross-section is symmetric with respect to the 
z-axis, the axial strain at^ = 0 and curvature about the z-axis are given by Equations 3.13. 
Substitution of Equations 3.46 and 3.47 into Equations 3.13 gives: 
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Hence, beam curvature with no axial strain at = 0 can be obtained by setting Vi = - V 2 
and V 3 = -V4. Similarly, pure axial strain with no beam curvature is achieved when Vt = V 2 
and V 3 = V 4 . 



C. ACTUATOR EFFECTIVENESS 

The results Euler-Bemoulli models of actuator produced curvature for the single, 
symmetric pair, single stacked and symmetric stacked pair configurations were compared 
for different thicknesses of an aluminum beam. All actuator configurations were assumed 
to cover the entire width of the beam and the bond between the piezoceramics and the 
beam was assumed to be negligibly thin. Material properties for Navy Type II PZT were 
used for the piezoceramic and a maximum applied electric field of 20 volts/mil, a typical 
operating limit, was applied to all actuators. Figure 11 shows the curvature of the 
composite actuator/beam as a function of the thickness of the individual piezoceramic 
elements for a 0.79 mm (0.031 in) thick beam. It can be clearly seen from Figure 1 1 that 
the symmetric actuator configurations are more effective at producing beam curvature 
than those with actuators on one side of the beam. This is due to the concentration of the 
piezoceramic at the outer edges of the composite beam for the symmetric configurations, 
where the actuator stress contributes most effectively to a bending moment about the 
central axis of the beam. An optimum piezoceramic thickness for each actuator type is 
also evident. Thicker than optimum actuators contribute more to the stiffhess of the 
composite beam cross-section than to actuator bending moment and hence result in lower 
than optimum curvature. The optimum thicknesses for the two stacked configurations 
here are exactly half those of the corresponding single actuator configurations due to the 
assumption of negligible bond thickness. The stacked configurations therefore require half 
the input voltage of the single actuator configurations to achieve a given electric field 
across the piezoceramic to produce a desired level of actuation strain. 
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Figure 11. Beam Curvature vs. Piezoceramic Thickness and Actuator Type for 0. 79 mm 
Aluminum Beam and 20 Volts/mil Electric Field applied to Actuators 

Figure 12 depicts the performance of the same actuators as in Figure 11 for an 
aluminum beam with twice the thickness. The optimum thickness of the actuators for each 
configuration is found to correspondingly double, however, the attainable curvature of the 
1 .58 mm beam is only half that attainable for the 0.79 mm beam at the 20 volts/mil electric 
field limitation. This is due to the maximum actuation strain of the piezoceramic 
remaining constant. The limit on actuation strain for a given piezoceramic results in an 
inverse relationship between beam thickness and the curvature attainable using an actuator 
of optimum thickness. 

Figure 13 and Figure 14 show the same information as Figure 11 and Figure 12, 
respectively, but with the actuator thickness plotted as a percentage of the thickness of the 
aluminum beam. Figure 15 and Figure 16 show the effect of actuator width for a constant 
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actuator thickness of 0.26 mm and beam thicknesses of 0.79 mm and 1.58 mm 
respectively. Full width actuators are clearly the most effective in producing beam 
curvature for all actuator configurations. The reason for this can be seen by examining the 
curvature equation in Equations 3.26, or the equivalent equation for other actuator 
configurations. Actuator width appears as a factor of the entire numerator but of only one 
term in the denominator, given by Equation 3.23. Increasing actuator width therefore 
always results in an increase in beam curvature, provided all other terms remain constant. 
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Figure 12. Beam Curvature vs. Piezoceramic Thickness and Actuator Type for 1.58 mm 
Aluminum Beam and 20 Volts/mi I Electric Field applied to Actuators 
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Examination of Equation 3.26 also reveals the material properties desirable in a 
piezoceramic actuator. High material modulus Ep, piezoelectric coefficient dsi and 
operating voltage, or equivalently operating electric field, all result in increased actuator 
effectiveness. In general, the higher the product of these three properties, the more 
effective the actuator. 
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Figure 13. Beam Curvature vs. Piezoceramic Thickness as a Percentage of Beam 
Thickness and Actuator Type for 0. 79 mm Aluminum Beam and 20 Volts/mil Electric 

Field applied to Actuators 
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Figure 14, Beam Curvature vs. Piezoceramic Thickness as a Percentage of Beam 
Thickness and Actuator Type for L58 mm Aluminum Beam and 20 Volts/mil Electric 

Field applied to Actuators 
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Figure 15. Beam Curvature vs. Piezoceramic Actuator Width and Type for 0. 79 mm 
Aluminum Beam and 20 Volts/mil Electric Field applied to Actuators 
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Figure 16. Beam Curvature vs. Piezoceramic Actuator Width and Type for 1.58 mm 
Aluminum Beam and 20 Volts/mi I Electric Field applied to Actuators 
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IV. ACTUATOR PLACEMENT 



A. DEFORMATION OF EULER-BERNOULLI BEAM UNDER ACTUATOR 

BENDING MOMENT 



1. Transverse Deflection of a Beam Segment due to an Actuator Moment 

A section of laminated beam with uniform cross-sectional properties throughout its 
length will deform with uniform curvature over its length when subjected to a pure 
bending moment. For small displacements, defined as beam slope much less than 1, this 
behavior can be approximated by; 



d\ _ M, 

dx‘ (£7). 




<x£x, 



+ /. 



k 7 



dx 




(4,1) 



where x is the coordinate along the length of the beam, Xk is the x coordinate of the start of 
beam segment k, 4 is the x length of the beam segment, Vo(x) is the transverse 
displacement of the beam's neutral axis in the >^-direction, Mk is the bending moment 
applied to beam segment k about its neutral axis and 



(£ 7 ). = ZE,I. 



(4,2) 



is the flexural rigidity of the composite beam segment about its neutral axis. Et and /, are 
the Young's modulus and moment of inertia about the beam segment neutral axis of the /th 
lamina on segment k. 

The transverse displacement at any point x on beam segment k can be obtained by 
integrating Equation 4.1 twice over the length of the segment. Integrating once gives the 
slope of the beam segment as a function of x: 
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(x, £l£x.+/,) (4.3) 






Integrating a second time gives; 



M, 



^oW = ^o(^*) + - x,) + ~ ^ ^ ^ 0 (4-4) 



Thus, for small displacements, the transverse displacement along the length of a beam 
segment is a function of the displacement and slope at the start of the beam segment, the 
curvature of the beam segment due to the actuator moment and the position x-Xk on the 
beam segment. 



2. Transverse Deflection of Complete Beam due to Actuator Moments 




Figure 17. Beam with Attached Piezoceramic Actuators 

Figure 17 shows a cantilever beam with n piezoceramic actuators bonded to the 
beam at locations x, (/' = 1,2,...,«). If the actuators are perfectly bonded to the beam over 
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their entire lengths and an actuator's bending moment is proportional to its input voltage, 
then each actuator will produce a beam curvature per unit voltage: 




(4.5) 



over its respective beam segment, where Ft is the input voltage of actuator k and k* is the 
curvature of the corresponding composite actuator^eam section. A suitable value for Kk 
for a given actuator/beam configuration can be obtained as discussed in Chapter II. 

If the deformation of the beam is static and no external loads are applied, then the 
deformation of the beam is strictly a function of the actuator locations, lengths, 
actuator^eam section curvatures per unit voltage and actuator input voltages. The static 
deformation of the beam on actuator k can be obtained by substituting Equation 4.5 into 
Equations 4.3 and 4.4, giving: 



Similarly, since the sections of the beam between the actuators have no applied 
moments, and hence zero curvature, the static deformation of the beam between actuator 
k and k-^1 can be obtained fi'om Equations 4.3 and 4.4 with the moment terms set to zero. 
Thus, 



Vo'U) = Vo'(jCt) + ^*Ft(x-x,) (x, <x<x, +/,) 



(4.6) 



V. 



-o(x) = v„(x,) + Vo’(x,X^- 




(x-xj' (x, <x<x*+/J (4.7) 



Vo'W = ^o'(^t+4) 



(4.8) 



^0 W = + 4 ) + V, ' (x, + /, )[x - (x^ + /^)] (x, + /, < X < x^,, ) (4.9) 
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Equations 4.6 through 4.9 include the terms vo(xk), vo(xk), vd(xk-^lii) and vo'(xk-^lk), 
which depend on the effects of actuators 1 through k-1 but not actuators k+ 1 through ii. 
Thus, a complete expression for the transverse deflection of the beam as a function of x 
can be obtained using Equations 4.6 through 4.9 by starting at the root of the beam (x=0) 
and working outward to the tip (x=L). 

Replacing the variable Vo with y, the slope y’ and displacement y of the beam are 
both zero prior to the first actuator; 

> (O < X < X,) 

The slope and displacement of the beam for x on the first actuator are functions of 
X and the properties of the first actuator only: 

>'’(x) = >^'(x,) + /:,F,(x - X,) 

= 0 + /CjF,(x- X,) 

■=Kr,(x-x,) 

yix)=y{x,) + y'(x,){x-x,) + ^^(x-xy' 

K.V. / V 
= 0 + 0 + -^(x-x,) 



The slope for x between the first and second actuators is the slope achieved at the 
end of the first actuator while the displacement for x in this range is the displacement at 
the end of the first actuator plus the slope at the end of the first actuator times the distance 
from the end of the first actuator; 



(x, < x< X, +/,) 



yw = ^ = o 

ax 

y(x) = 0 
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y(x) = y(x,+/,) 

>'(x) = j(x, +/,) + y(x, +/,)[x-(x, +/,)] 

= ^// + A:y./.[x-(x.+/.)] 



^ 2 



(x, +/, <X<xJ 



Continuing in this manner leads to the following general slope and 
displacement relations for a beam with n attached piezoceramic actuators; 

a) For X prior to the first actuator 



yw = ol 

, \ (0<x<x.) 

t(x)=OJ ' 



b) For X on the first actuator 



y(x) = ^,y(x-x,) ' 

^(x) = ^^(x-x,y 



(x, < X < X, + /,) 



c) For X on actuator k 



y(x) = ^*F*(x-xJ+S^y/, 

y{x) = -x,Y + Z [x-x,- 



X* <x<x, +4 
U = 2,3,...,n , 



d) For X between actuator k and k+1 



y(x)=Y.K.v.i. 

7=1 



Xk 

y = l,2,...,n-l , 



transverse 



(4.10) 



(4.11) 



(4.12) 



(4.13) 
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e) For X between actuator n and the free end of the beam (x = lu) 



n 










(4.14) 



B. OPTIMIZATION PROBLEM IN KUHN-TUCKER FORM 

It is desired to determine the optimum placement of piezoceramic actuators to 
approximate a desired beam shape function under actuator-induced deformation. 
Optimization requires a suitable cost function to serve as a measure of performance for a 
given actuator configuration. The optimum actuator configuration achieves the absolute 
minimum value of the cost function among all possible values of the input variables. For 
piezoceramic actuators, these input variables can be expressed as actuator position, length, 
input voltage, and beam curvature per unit actuator input voltage. 

An appropriate cost function for measuring actuator performance in static 
deflection of a cantilever beam is: 



where x is the coordinate along the length of the beam and y(x) and yd(x) are the actual 
and desired transverse beam deflections at each point along the beam, respectively, as 
shown in Figure 17. Values of x range from 0 at the constrained end to L at the 



variables that produces a y(x) which provides a continuous "least squares" fit to the 
desired shape function over the length of the beam. 

The actuator locations and lengths are constrained by; 



L 




(4.15) 



0 



unconstrained end of the beam. Equation 4.15 is minimized by the combination of input 
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(4.16) 



X, > 0 

x„ + l„^L 



l,>0 {i=\,2,. ,n) 



where n is the number of actuators and x, and // are the locations and lengths of the 



actuators as shown in Figure 17. Equations 4.16 ensure that all actuators are located 
entirely on the beam with no overlap and that all have nonnegative length. Constraints on 
the maximum actuator lengths can also be included. Furthermore, the input voltages to 
the piezoceramic actuators must lie within an operating range specified by; 



where and V„ax are the actuator's minimum and maximum operating voltage and ¥ is 
the input voltage to the rth actuator. 

The problem of optimizing Equation 1 within the constraints of Equations 4.16 and 
4.17 can be expressed in Kuhn-Tucker (KT) form [Ref 8] as; Minimize Equation 4.15 
subject to: 




(4.17) 



gi(z) = -x, ^0 

g;(z) = Xi_, -Xj <0 (i = 2,3,...,n) 

= +1„ -L<0 

gn.ui(z) = -li<0 (i = l,2,...,n) 

g 2 „.ui (z) = (X - \W, - V_ ) < 0 (i = 1,2, . . . , n) 



(4.18) 



where the vector of input variables, z, is given by 




(4.19) 
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A KT problem is typically solved by constructing a Lagrangian, given by; 



L = J(z) + //,g,(z) + 



(4.20) 



and solving the system of equations: 



dz 



Mi) 






^ik) , ^zk) . . ^mk) 



= Q 



1=0 ifg,k)<0 
[>0 ifg.(z) = 0 



(/ = u 




(4.21) 



where m is the number of constraint equations. Equations 4.21 with all /i, terms set to 0 
may have no solutions or one or more solutions for the values of the input variables, with 
each solution corresponding to a local minimum or maximum of Equation 4.15. A local 
minimum may or may not be the global minimum, if a global minimum exists. If Equation 
4.15 has either no finite global minimum or a global minimum outside the constrained 
range of input variables, the minimum value of Equation 1 over the allowable range of 
input variables may lie on a constraint boundary. The /i, terms in Equations 4.20 and 4.21 
allow for the possibility that the minimum value of the cost function in Equation 4.15 over 
the bounded range of possible input variables specified by Equation 4.18 lies on one or 
more constraint boundaries. The minimum value of the cost function on constraint 
boundary / can be found by including //, as an additional variable to be solved for in the z 
vector and adding the corresponding constraint equation gi=0 to the system of equations 
to be solved in Equations 4.21. Differentiating Equation 4.20 with respect to the 
additional variable /// can be seen to return the appropriate constraint equation. 
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C. SELECTED SHAPE FUNCTIONS FOR ANALYSIS 

An analytical expression for the cost function to be minimized can be obtained by 
substituting Equations 4.10 through 4.14 into Equation 4.15, provided that all terms 
involving the desired beam shape function yd(x) are integrable. Two suitable shape 
functions were selected for evaluation of optimal actuator placement. The first, a 
parabolic shape function given by: 

y,{x) = Cx^ (4.22) 

was selected as a due to its mathematical simplicity and its applicability to common 
structures such as antennas. This shape function has uniform curvature C over the length 
of the beam. A second shape function given by: 

y^ix) = C[l - cos(wx)] (4.23) 

where m = n/L, was chosen to analyze the effect of non-uniform desired curvature and an 
inflection point on the optimization process. This function has curvature: 

y<i = 

which is positive over the first half of the beam and negative over the second half The 
inflection point for this function is at x = 0.5L, where the curvature is zero. 

Substituting Equation 4.22 and Equations 4.10 through 4.14 into Equation 4.15 
gives an expression for the cost function associated with the parabolic desired shape: 
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5 






C^x, 



\^2t w C ^/^U A 



+ S^ + 
*=1 



U/'_ _ 2l(^*+^*)'-V 

[2(^2 it ^)^0k ~^^\k J 2 

+lkf -Xj 



^^Xk^Ok - ~^~^A*. 



0^ A: 



n-\ 



c 



2 ^*^1 



— ^ -2C6„ 



+ E^(-2C*.. 



k=\ 



■^^Ot (^*+1 ^k ^i) 



l’-(x +/) r-(x+/„) 

+ C2 — —-2Cb, ^ 

5 >" 4 

+ (- 2C6o„ + 6,/ ) ^3 + 26. A 

+ bj(L-x„-l„) 



'\n^0n 



(4.25) 



where the coefficients ay and by are obtained by rewriting the beam displacements in 
Equations 4. 10 through 4. 14 as polynomials in x; 



For X prior to the first actuator; 

>'(x) = 0 (o<x<x,) 



(4.26) 



For X on actuator k. 



j(x) = a2*x" +a„x + ao, 



"x, <x<x, +/, 

Kk = 1,2,...,« 



(4.27) 
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For X between actuator k and k+1: 



+/* 

Kk = \,l,...,n-\ j 



Ax) = b,,x + b^, 



(4.28) 



For X between actuator n and the free end of the beam (x = L) : 

yix)^b,„x + b^„ (x„+/„<x<l) (4.29) 



giving: 



^2k ~ 






(A = 1,2, ...,«) 






(* = 0 

-K,V,x,^^Kyj, {k = l,3,...,n) 



-K,V,x, 

k-\ 



«01 



2 

2 k-\ 



(A: = l) 



K,V,x,^ ViTF/f J~\ a \ 

2 + 2 ) = 2,3,...,w) 



*u=ZW. {k = l,2,...,n) 



bok =-i:Kr,h 



1=1 



X, + 

^ ^ 2, 



(* = l^,...,n) 



(4.30) 



An equation for the cost function associated with the second shape function can 
similarly be obtained by substituting Equation 4.23 and 4.26 through 4.30 into Equation 
4.15, giving: 
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J = c 



X, - — sin(/wx, ) + — + — sin(2wL) 
' /w ^ 2 4/w ^ ’ 



+E 

Jt=l 



a 






Ik 






i^k +hy -X,' 



2] (^k+ikf -X, 



+ 2 «u(«ot 









^^2kC 



+ 2a2f^C 



[(^* +^t)cos(/n(x, +/*))-x, cos(/nx,)] 

sin(/w(x, + ^* )) - f ^ sin(mx, 



{^k + ) 2 



^ /W /W y 



+ ^^[cos(/w(x, +/*))- cos(/wx, )] 

^[(^* + Ik ) sin(/n(x, +/*))- x* sin(/wx* )] 
^ l ^ |sin(w(x, +/,))-sin(wx,)] 



m 
2 ^C 
m 



n-\ 

*=1 



2 ^k^X 



'\k 






+ (^* -C)'(x*,, -X* -/,) 



2^uC* 

26 „C 



[cos(/nx,^,)-cos(/w(x, +/,))] 



^[x,„ sin(/wx*^,)-(x, +/,)sin(/w(x, +/* 
[sin(mx,^i)- sin(w(x, + /, ))] 



(continued on next page) 
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+ b,: ^ * 2b„(K. - c)— 

+ [cos(wi) - cos(m(>:. + /.))] 

2jy ^ 

+ —^\Lsm{mL)-[x„ +/„)sin(/w(x„ +/„) 

- sin(w(x„ + /„))] 



(4.31) 



D ALGORITHM FOR SOLUTION OF OPTIMIZATION PROBLEM FOR 
ACTUATORS WITH FIXED DIMENSIONS AND PROPERTIES 

The MATLAB software package was used for development of an algorithm to 
optimize actuator placement and voltage for a given shape function and fixed actuator and 
beam dimensions and properties. The shape functions of Equations 4.22 and 4.23, with 
Equations 4.25 and 4.31, respectively as corresponding cost functions, were used for 
evaluation. Inspection of Equation 4.25, with appropriate substitutions from Equations 
4.30, reveals that its terms are up to fifth order polynomials in the actuator locations x, and 
lengths /, multiplied by second order polynomials in the actuator curvatures 
Equation 4.31 is similar, with sine and cosine terms replacing some of the polynomial 
terms. If the actuator dimensions and properties are fixed, then the lengths /, and 
curvatures per unit input voltage will be constant values. Either cost function is therefore 
minimized by optimizing the remaining variables x, and Vi. Equivalently, if only the 
actuator lengths are specified, the cost function is minimized by optimizing the actuator 
locations x, and curvatures xj = KiVi. The actuator configuration and cross-sectional 
properties can then be chosen to maximize the actuator curvature per unit voltage as 
discussed in Chapter III. 

The MATLAB Optimization Toolbox includes three functions intended for the 
optimization of multivariable functions. Two functions, fminu and jmins, perform 
unconstrained optimization of the input variables of a cost function, while one function. 
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constr, performs optimization of the input variables subject to constraints on their values. 
Fmins uses the simplex search algorithm developed by Nelder and Mead. Fminu uses the 
quasi-Newton method of Broyden, Fletcher, Goldfarb and Sharmo (BFGS) with a mixed 
quadratic and cubic line search procedure to determine a search direction at each iteration. 
Constr also uses a BFGS quasi-Newton method to solve a quadratic subproblem at each 
iteration of a sequential quadratic programming routine. Grace [Ref 9] provides a 
thorough description of the BFGS and sequential quadratic programming algorithms but 
omits detailed discussion of the Nelder and Mead technique. 

Initial attempts to optimize actuator placement and voltages simultaneously using 
constr and fminu proved unsuccessfiil. These algorithms returned unreliable results with 
minimal change in the initial actuator locations and voltages specified for each test, 
regardless of the settings of the tolerances on the input variables and cost function 
evaluations and other options associated with the functions. Fmins produced much better 
results in general but was found to return erroneous results for some sets of initial 
conditions for three or more actuators. The likelihood of erroneous results increased 
directly with the number of actuators evaluated. A more robust method of optimizing 
actuator locations was therefore judged to be needed. 

A revised approach to the optimization problem was developed, drawing on the 
work of Clark and Fuller [Ref 10] and Wang, Burdisso and Fuller [Ref 11] to optimize 
actuator locations for structural acoustic control. The latter group noted that the 
mathematical inconsistency in order of the voltages and actuator locations in their 
objective function formulation requires independent solution for the voltages and 
locations. Both groups therefore determined optimum actuator voltages for each set of 
actuator locations evaluated in their algorithms to determine optimum actuator placement. 
A similar approach can be applied to this problem by using an outer optimization 
algorithm to determine the actuator locations which minimize Equation 4.25 or 4.3 1 with 
an embedded algorithm to determine the optimum actuator input voltages for each set of 
actuator positions evaluated. The cost function for the outer algorithm thus becomes a 
fifth order polynomial in the actuator locations x,, with sine and cosine terms if Equation 
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4.31 is used, while the cost function for the embedded voltage optimization algorithm 
becomes a second order polynomial in the actuator input voltages F,. 

A two-stage solution algorithm to determine the optimum actuator locations and 
voltages for a given set of actuator and beam dimensions and properties was implemented 
as described above using MATLAB. The algorithm was tested using either constr,fminu 
or fmins to perform the required optimization. Both comtr and fminu were again found to 
return unreliable results with almost no change in the specified initial conditions but the 
fmins simplex search routine was found to be highly reliable, converging to the same 
results for a given number of actuators and set of actuator and beam dimensions and 
properties, over a wide range of specified initial actuator positions. This is consistent with 
the observations of Parkinson and Hutchinson [Refs. 12, 13], who noted that the Nelder 
and Mead Simplex (NMS) algorithm has proven to be robust, although less efficient than 
some unconstrained optimization algorithms. A block diagram of the algorithm 
incorporating ^/>25 is shown in Figure 18. The outer stage of the algorithm uses finins to 
perform a simplex search to determine the actuator positions, starting from an initial 
estimate, which minimize the cost function to a set tolerance. The inner stage of the 
algorithm uses fmins to determine the optimum actuator input voltages for each set of 
actuator positions evaluated by the outer stage. The effect of actuator input voltage is 
therefore completely transparent to the outer stage, while the actuator locations are fixed 
for each evaluation of the inner stage. 

The use of the unconstrained optimization algorithm fmins to solve the constrained 
actuator location and voltage optimization problems is acceptable provided that the 
optimum solution is at the global minimum of the cost function and does not violate any of 
the constraints. Since the integrand of Equation 4.15 is always non-negative, the value of 
the cost function will also be non-negative and the cost function will have a global 
minimum. The algorithm presented here is therefore suitable for shape functions which do 
not result in actuator overlap and which require actuator input voltages within the 
allowable range. Complex desired shape functions with significant curvature variation 
over distances on the order of an actuator length or less might result in solutions with 
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actuator overlap, requiring that the NMS routine be modified to include constraints or 
replaced with another suitable constrained optimization algorithm for both voltage and 
position optimization. 




Figure 18. Actuator Placement Optimization Algorithm 



E. VALroATION OF OPTIMIZATION ALGORITHM 
1. Results for Parabolic Shape Function 

The optimization algorithm was implemented using the MATLAB program 
actplace.m and its associated subroutines, which are presented in the appendix. The 
program was tasked with determining the optimum locations for between one and five 
actuators on a one meter long cantilever beam to best approximate a desired parabolic 
beam shape with a tip displacement equal to 0.1% of the beam length. The beam was 
specified to have the properties of 7075-T6 aluminum, a length of 1 meter, a width of 5.08 
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cm (2.00 in.) and a thickness of 1.575 mm (0.062 in ). The beam curvature per unit 
actuator input voltage was calculated for a symmetrically mounted pair of 0.26 mm (10 
mil) thick Navy Type II piezoceramic actuators covering 75% of the width of the beam. 
Each actuator length was fixed at 10% of the overall beam length. The algorithm was 
tested wdth at least four sets of arbitrarily selected initial actuator locations for each 
number of actuators to ensure convergence to the same optimum values. The initial 
configurations tested included all actuators clustered at the root of the beam, all actuators 
clustered at the center of the beam and actuators distributed over the length of the beam. 
An additional test was performed with the initial actuator locations within a few percent of 
beam length of the optimum locations found on previous trials for comparison of the 
number of floating point operations required to reach convergence. Initial actuator 
voltages were the optimum voltages for the initial actuator locations for each case, 
determined usmgjmim, except where noted. 

The results for each number of actuators are tabulated with respect to each set of 
initial conditions in Table 2 through Table 6. The optimization algorithm proved to have 
minimal sensitivity to the initial actuator locations specified for up to five actuators, with 
two exceptions denoted by asterisks in Table 5 and Table 6. Grouping the initial locations 
for four actuators adjacent to one another at the center of the beam and using optimum 
voltages for this initial configuration as an initial estimate for all voltage evaluations 
prevented the inner voltage optimization algorithm fi-om converging within the set 
maximum of 10,000 iterations at some point during the optimization process. Using the 
same actuator locations with initial voltage estimates of 0 produced erroneous results with 
some actuator overlap, as shown in the Run 2 column of Table 5. Similarly, grouping the 
initial locations for five actuators adjacent to one another at the root of the beam produced 
erroneous results with the first actuator overhanging the fixed end of the beam, as shown 
in the Run 1 column of Table 6. For both of these cases, slightly changing the initial 
conditions to allow 1 to 2% of beam length spacing between the initial actuator positions 
produced results consistent with those for the remaining sets of initial conditions tested, as 
shown in the Run 5 columns of the respective tables. The optimum actuator positions 
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obtained in the remaining cases for each number of actuators differed by at most lO"^, 
which was the set tolerance on the actuator locations for termination of the fmim function. 
The optimum deformed beam shape is plotted versus the desired parabolic shape for each 
number of actuators in Figure 19 through Figure 23. 

The results of the actuator voltage optimization algorithm were compared to the 
results of exhaustive searches of all possible combinations of discrete sets of actuator 
voltages, for the one, two and three actuator cases. The actuator voltages were varied 
discretely over specified ranges for actuator locations fixed at the optimum locations 
found by the placement optimization algorithm. Figure 24 through Figure 27 show plots 
of the variation of error with actuator input voltage. For all cases, the error due to the 
voltages obtained from the optimization was less than the smallest error obtained by the 
discrete search. Additionally, the voltages obtained from the optimization algorithm 
coincide with the locations of the minimum error in the corresponding plots. 

The results of the actuator placement optimization algorithms were similarly 
compared to the results of exhaustive searches of all possible combinations of discrete sets 
of actuator locations for the one, two and three actuator cases. The actuator locations 
were varied discretely over specified ranges and the optimum actuator input voltages were 
determined for each combination of actuator locations using the voltage optimization 
algorithm. Figure 28 through Figure 3 1 show plots of the variation of error with actuator 
location. As for the voltages, the error due to the actuator locations obtained from the 
optimization was less than the smallest error obtained by the discrete search and the 
locations obtained from the optimization algorithm coincide with the locations of the 
minimum error in the corresponding plots. Table 7 shows the results for three actuators 
which can be compared to the results of the optimal placement algorithm in Table 4. 

The comparison of the slopes of Figure 25 through Figure 27 and Figure 29 
through Figure 31 illustrates the relative importance of each of the actuators in 
determining the error between the actual and desired beam shapes. The location and input 
voltage, or equivalently curvature, of actuator 1, the actuator closest to the cantilevered 
end of the beam have the most significant effect. Similarly, the location and input voltage 
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of actuator 2 have a greater impact than those of actuator 3. This is consistent with the 
fact that each actuator affects the shape of the beam from its location to the beam’s free 
end but has no effect on the shape of the beam between the cantilevered end and its 
location. Actuator 1 therefore affects the greatest portion of the beam’s length, followed 
by actuator 2, then 3, etc. 

2. Results for Second Shape Function 

Following evaluation of the optimization algorithm’s performance with the 
parabolic shape function, the algorithm was tasked with determining the optimum location 
for between one and four actuators to best approximate the beam shape given by Equation 
4.23 with a beam tip displacement again equal to 0.1% of the beam length. The beam and 
actuators were given the same dimensions and properties as in the previous tests. The 
algorithm was tested with at least four sets of initial actuator locations for each number of 
actuators to check convergence to the same optimum values. The initial configurations 
tested again included all actuators clustered at the root of the beam, all actuators clustered 
at the center of the beam, actuators distributed over the length of the beam and initial 
locations within a few percent of beam length of the optimum locations found on previous 
trials. Initial actuator voltages were the optimum voltages for the initial actuator locations 
for each case. 

The results of the tests for the shape function of Equation 4.23 are presented in 
Figure 32 through Figure 36 and Table 8 through Table 11. The optimization algorithm 
showed greater sensitivity to initial conditions for this shape function than for the 
parabola, particularly for the four actuator test, but again converged to essentially the 
same results for a variety of sets of initial conditions for each number of actuators. 
Several sets of initial conditions prevented the inner voltage optimization algorithm from 
converging within the set maximum of 10,000 iterations at some point during the 
optimization process. Grouping the initial actuator positions at the root of the beam 
resulted in the algorithm converging on local minima of the cost function for the two and 
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four actuator cases, as shown in the Run 1 column of the respective tables. While the 
result of the two actuator case for these initial conditions, with the first actuator 
overhanging the fixed end of the beam, is obviously erroneous, the result of the four 
actuator case is a valid local minimum. Comparison of the actuator locations which result 
in the global minimum of the cost function in Figure 35 and the actuator locations 
resulting in the local minimum in Figure 36 illustrates this. 

The curvature of the desired shape function, given by Equation 4.24, is 
decreasingly positive from the root to the midpoint of the beam and then increasingly 
negative from the midpoint to the free end of the beam. It can thus be qualitatively seen 
that actuators must be located in the regions of high desired curvature near the ends of the 
beam to most closely approximate the desired shape. Similarly, actuators in the low 
curvature region near the midpoint of the beam provide relatively little benefit in reducing 
the cost function. The combination of three actuators producing successively lower 
curvatures near the root and one actuator producing a high curvature in the opposite 
direction near the free end of the beam in Figure 36 provides a relatively good 
approximation of the shape function over its high-curvature portion near the root of the 
beam but a poorer approximation of the shape function over its high-curvature portion 
near the free end. In contrast, the combination of two actuators near the root and two 
actuators near the free end of the beam in Figure 35 provides a slightly poorer 
approximation of the shape function over its high-curvature portion near the root of the 
beam but a much better approximation of the shape function over its high-curvature 
portion near the free end, resulting in a lower overall value of the cost function than the 
arrangement in Figure 36. Grouping the initial actuator positions near the root of the 
beam prevented the optimization algorithm from moving more than one actuator past the 
inflection point in the desired shape function at x = 0.5£. However, the initial conditions 
for runs 4 through 6 in Table 1 1, which start two actuators on each side of the inflection 
point, allow the optimization algorithm to achieve the correct optimum positions. 
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3 . 



General Observations 



The two shape functions considered here illustrate the importance of properly 
selecting initial conditions for the optimization algorithm. Allowing some space between 
initial actuator locations was proven to be beneficial, but not always required to obtain 
convergence. The presence of inflection points or significant variations in the curvature of 
the desired shape function can result in actuator locations producing a local minimum in 
the error cost function and must therefore be considered in selecting initial conditions. 
Knowledge of the approximate locations of the optimum actuator positions should be used 
if possible, and multiple sets of initial conditions should be evaluated to ensure the 
minimum obtained is the global minimum. Exhaustive searches of discrete sets of actuator 
locations can be helpful in selecting initial conditions, but such searches become 
increasingly computationally intensive as the number of actuators is increased. 
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Actual v«. Doslrod Beam Shapes 




Actual 

Desired 



Figure 19. Optimal Location for One 10% Length Actuator for Parabolic Desired Shape 



Run Number 


1 


2 


3 


4 


Initial Actuator 
Location (m) 


0 


0.2 


0.6 


0.25 


Initial Actuator 
Voltage (V) 


25.23 


36.71 


104.59 


40.74 


Optimum Actuator 
Location (m) 


0.2419 


0.2420 


0.2419 


0.2419 


Optimum Actuator 
Voltage (V) 


40.05 


40.05 


40.05 


40.05 


Error (10'^) 


1.335 


1.335 


1.335 


1.335 


Operations (10^) 


3.08 


1.53 


2.93 


1.10 



Table 2. Optimization Algorithm Results for One 10% Length Actuator for Parabolic 

Desired Shape 
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X 10 



Actual vs. Ddttlrdd Beam Shapes 





Figure 20. Optimal Locations for Two 10% Length Actuators for Parabolic Desired 

Shape 



Run Number 


1 


2 


3 


4 


Initial Actuator 


0.0 


0.4 


0.25 


0.1 


Locations (m) 


0.1 


0.5 


0.65 


0.5 


Initial Actuator 


-26.64 


94.79 


37.66 


20.85 


Voltages (V) 


61.55 


-50.78 


12.97 


26.95 


Optimum Actuator 


0.1220 


0.1220 


0.1220 


0.1220 


Locations (m) 


0.5394 


0.5394 


0.5394 


0.5394 


Optimum Actuator 


23.50 


23.50 


23.50 


23.50 


Voltages (V) 


24.85 


24.85 


24.85 


24.85 


Error (10'*“) 


1.366 


1.366 


1.366 


1.366 


Operations (10^) 


2.88 


1.72 


1.55 


0.83 



Table 3. Optimization Algorithm Results for Two 10% Length Actuators for Parabolic 

Desired Shape 
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Beam DleplacementY (m) 



Actual v«. Daclred Beam Shapes 




Actual 

Desired 



X(m) 

Figure 21. Optimal Locations for Three lOVo Length Actuators for Parabolic Desired 

Shape 
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Run Number 


1 


2 


3 


4 


5 


6 


Initial Actuator 


0.0 


0.35 


0.20 


0.1 


0.05 


0.075 


Locations (m) 


0.1 


0.45 


0.45 


0.4 


0.35 


0.370 




0.2 


0.55 


0.70 


0.7 


0.65 


0.665 


Initial Actuator 


16.61 


100.09 


31.45 


19.40 


14.69 


16.89 


Voltages (V) 


-40.88 


-99.82 


4.58 


17.62 


19.48 


18.18 




62.52 


50.20 


19.22 


17.94 


18.51 


18.25 


Optimum 


0.0728 


0.0728 


0.0728 


0.0728 


0.0728 


0.0728 


Actuator 


0.3678 


0.3678 


0.3678 


0.3678 


0.3678 


0.3678 


Locations (m) 


0.6641 


0.6641 


0.6642 


0.6641 


0.6642 


0.6641 


Optimum 


16.67 


16.67 


16.67 


16.67 


16.67 


16.67 


Actuator 


18.32 


18.32 


18.32 


18.32 


18.32 


18.32 


Voltages (V) 


18.30 


18.30 


18.30 


18.30 


18.30 


18.30 


Error (1 O'") 


2.779 


2.779 


2.779 


2.779 


2.779 


2.779 


Operations 

(■o') 


20.3 


8.63 


8.39 


4.48 


5.48 


3.17 



Table 4. Optimization Algorithm Results for Three 10% Length Actuators for Parabolic 

Desired Shape 



Run Number 


1 


2* 


3 


4 


5 


Initial Actuator 


0.0 


0.3 


0.1 


0.05 


0.28 


Locations (m) 


0.1 


0.4 


0.3 


0.25 


0.39 




0.2 


0.5 


0.6 


0.50 


0.51 




0.3 


0.6 


0.8 


0.75 


0.62 


Initial Actuator 


-22.33 


0 


17.51 


12.58 


64.31 


Voltages (V) 


70.51 


0 


12.75 


13.23 


-44.41 




-47.54 


0 


17.92 


16.22 


15.55 




38.53 


0 


8.98 


14.35 


16.00 


Optimum Actuator 


0.0462 


0.0794 


0.0462 


0.0462 


0.0462 


Locations (m) 


0.2742 


0.3700 


0.2743 


0.2742 


0.2743 




0.5038 


0.6141 


0.5039 


0.5038 


0.5039 




0.7334 


0.7697 


0.7335 


0.7334 


0.7335 


Optimum Actuator 


12.95 


17.35 


12.95 


12.94 


12.95 


Voltages (V) 


14.19 


16.87 


14.20 


14.19 


14.20 




14.20 


12.17 


14.19 


14.20 


14.20 




14.18 


9.74 


14.18 


14.18 


14.18 


Error (10-*^) 


7.528 


16.48 


7.528 


7.528 


7.528 


Operations (10^) 


7.88 


7.83 


2.69 


1.78 


5.27 



Table 5. Optimization Algorithm Results for Four 10% Length Actuators for Parabolic 

Desired Shape 



65 



B«am Dbplaeainonl Y (m) 



Actual vs. Ddsirdd Baam Shapes 




Actual 

Odsirad 



X{m) 

Figure 22. Optimal Locations for Four 10% Length Actuators for Parabolic Desired 

Shape 
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B«am Dkplaeomont Y (m) 




X(m) 

Figure 23. Optimal Locations for Five 10% Length Actuators for Parabolic Desired 

Shape 
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Run Number 


1* 


2 


3 


4 


5 


Initial Actuator 


0.0 


0.25 


0.05 


0.05 


0.00 


Locations (m) 


0.1 


0.35 


0.25 


0.20 


0.11 




0.2 


0.45 


0.45 


0.40 


0.22 




0.3 


0.55 


0.65 


0.60 


0.33 




0.4 


0.65 


0.85 


0.80 


0.44 


Initial Actuator 


-1.54 


27.22 


11.59 


11.25 


11.33 


Voltages (V) 


31.37 


12.93 


16.55 


6.56 


-2.36 




-30.73 


13.58 


4.06 


25.93 


3.04 




26.74 


-16.96 


23.07 


-3.48 


17.09 




17.24 


8.83 


-7.16 


12.04 


14.88 


Optimum Actuator 


-0.0100 


0.0298 


0.0298 


0.0298 


0.0298 


Locations (m) 


0.1531 


0.2155 


0.2155 


0.2155 


0.2156 




0.3548 


0.4031 


0.4031 


0.4031 


0.4031 




0.5575 


0.5907 


0.5906 


0.5907 


0.5906 




0.7610 


0.7782 


0.7782 


0.7782 


0.7782 


Optimum Actuator 


6.39 


10.62 


10.62 


10.62 


10.62 


Voltages (V) 


12.38 


11.60 


11.60 


11.60 


11.60 




12.51 


11.60 


11.60 


11.60 


11.59 




12.56 


11.60 


11.60 


11.59 


11.59 




12.57 


11.58 


11.58 


11.58 


11.58 


Error (10'*^) 


3.688 


2.310 


2.310 


2.310 


2.310 


Operations (10^) 


15.7 


14.5 


7.52 


6.80 


24.5 



Table 6. Optimization Algorithm Results for Five 10% Length Actuators for Parabolic 

Desired Shape 



68 



Error 



Error vs. Actuator Input Voltage for 1 Actuator 




Figure 24. Variation of Error with Input Voltage for One Actuator at Optimum Location 

for Parabolic Desired Shape 
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Error vs. Actuator Input Voltages for 2 Actuators 



J£ 10 



Run 1 

Actuator Locations; 0 .122 m Q.S394 m 




Figure 25. Variation of Error with Input Voltages for Two Actuators at Optimum 
Locations for Parabolic Desired Shape 
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^ Errsr vs. Actuaisr 1 a 2 Voltage atOptfinum Actuator 3 Voltago 




Figure 26. Variation of Error with Actuator 1 and 2 Input Voltages for Three Actuators 
at Optimum Locations and Optimum Actuator 3 Input Voltage for Parabolic Desired 

Shape 
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Error vs. Actuator 2 A 9 Voltage at OpUmum Actuator 1 Voltage 



jc IQ* 




Figure 27. Variation of Error with Actuator 2 and 3 Input Voltages for Three Actuators 
at Optimum Locations and Optimum Actuator 1 Input Voltage for Parabolic Desired 

Shape 
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Error 



Error vs. Actuator Location for 1 Actuator 




Figure 28. Variation of Error with Actuator Location for One Actuator for Parabolic 

Desired Shape 
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Error vs. Actuator Loeailons for 2 Actuators 




Figure 29. Variation of Error with Actuator Location for Two Actuators for Parabolic 

Desired Shape 



lA 



Run Number 


1 


2 


3 


4 


5 


6 


Actuator 


0.0- 


0.0- 


0.00- 


0.000- 


0.05- 


0.055- 


Location 


0.7 


0.4 


0.20 


0.100 


0.10 


0.095 


Ranges (m) 


0.1- 


0.2- 


0.25- 


0.325- 


0.35- 


0.355- 




0.8 


0.6 


0.45 


0.425 


0.40 


0.395 




0.2- 


0.5- 


0.55- 


0.625- 


0.65- 


0.655- 




0.9 


0.9 


0.75 


0.725 


0.70 


0.695 


Spacing 
Between 
Discrete 
Actuator 
Locations (m) 


0.1 


0.05 


0.025 


0.0125 


0.00625 


0.005 


Discrete 


0.1 


0.05 


0.075 


0.075 


0.07500 


0.075 


Actuator 


0.4 


0.35 


0.375 


0.375 


0.36875 


0.370 


Locations for 
Minimum Error 
(m) 


0.7 


0.65 


0.675 


0.675 


0.66625 


0.665 


Optimum 


19.40 


15.88 


16.96 


16.96 


16.87 


16.89 


Actuator 


17.62 


19.48 


18.68 


18.68 


18.06 


18.18 


Voltages (V) 


17.94 


18.50 


18.13 


18.13 


18.28 


18.25 


Error (x 10'") 


3.674 


3.210 


2.806 


2.806 


2.786 


2.783 



Table 7. Results of Discrete Actuator Placement Algorithm for Three Actuators for 

Parabolic Desired Shape 
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Erroi vs. Acaiator 1 a 2 LocaHon at Optimum Actuator 3 Location 




Figure 30. Variation of Error with Actuator 1 and 2 Location for Three Actuators with 
Actuator 3 at Optimum Location for Parabolic Desired Shape 



76 



Error vs. Actuator 2 a 3 Location atOptfcnum Actuator 1 Location 
Run 2 

Actuator Locations; Q.QSm Q.3S m Q.65in 

X IQ 

1.2 
1 - 
Q.0 ^ 

Q.6 
Q.4 
Q.2 
Q 

Q.6 

Q.S 

Actuator 2 (m) ° ^ “ ® Actuator 8 (m) 

Figure 31. Variation of Error with Actuator 2 and 3 Location for Three Actuators with 
Actuator 1 at Optimum Location for Parabolic Desired Shape 
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Figure 32. Optimal Location for One 10% Length Actuator for Equation 4.23 Desired 

Shape 



Run Number 


1 


2 


3 


4 


Initial Actuator 
Location (m) 


0 


0.2 


0.6 


0.1 


Initial Actuator 
Voltage (V) 


35.29 


49.93 


126.64 


41.74 


Optimum Actuator 
Location (m) 


0.0631 


0.0631 


0.0632 


0.0631 


Optimum Actuator 
Voltage (V) 


39.19 


39.19 


39.19 


39.19 


Error (10'^) 


1.243 


1.243 


1.243 


1.243 


Operations (10^) 


4.47 


4.22 


5.90 


2.86 



Table 8. Optimization Algorithm Results for One 10% Length Actuator for Equation 

4.23 Desired Shape 



78 



Actual vs. D««lfdd Boa in Shapos 




Actual 

Dosirod 



Figure 33. Optimal Locations for Two 10% Length Actuators for Equation 4.23 Desired 

Shape 



Run Number 


1* 


2 


3 


4 


5 


Initial Actuator 


0.0 


0.4 


0.25 


0.1 


0.1 


Locations (m) 


0.1 


0.5 


0.65 


0.25 


0.75 


Initial Actuator 


13.64 


214.07 


67.89 


50.72 


44.21 


Voltages (V) 


25.70 


-190.42 


-54.75 


-12.17 


-31.68 


Optimum Actuator 


-0.0100 


0.0942 


0.0942 


0.0941 


0.0942 


Locations (m) 


0.0989 


0.7441 


0.7441 


0.7441 


0.7441 


Optimum Actuator 


12.08 


43.70 


43.70 


43.70 


43.70 


Voltages (V) 


27.27 


-29.32 


-29.31 


-29.31 


-29.30 


Error (1 O'*”) 


12.25 


1.698 


1.698 


1.698 


1.698 


Operations (10*) 


2.28 


3.96 


2.96 


3.41 


1.81 



Table 9. Optimization Algorithm Results for Two 10% Length Actuators for Equation 

4.23 Desired Shape 
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(ui) uieoB 



X IQ 



Actual VC. Doclred Beam Shapes 





Figure 34. Optimal Locations for Three 10% Length Actuators for Equation 4.23 

Desired Shape 
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Run Number 


1 


2 


3* 


4 


5 


Initial Actuator 


0.0 


0.35 


0.20 


0.1 


0.05 


Locations (m) 


0.1 


0.45 


0.45 


0.4 


0.20 




0.2 


0.55 


0.70 


0.7 


0.70 


Initial Actuator 


21.44 


195.57 


66.51 


43.49 


27.29 


Voltages (V) 


11.32 


-183.43 


-26.36 


3.27 


18.43 




6.28 


18.63 


-19.35 


-28.43 


-26.67 


Optimum 


0.0268 


0.0268 


Did Not 


0.0268 


0.0267 


Actuator 


0.2133 


0.2133 


Converge 


0.2133 


0.2132 


Locations (m) 


0.7184 


0.7184 




0.7184 


0.7184 


Optimum 


24.60 


24.60 




24.61 


24.59 


Actuator 


21.18 


21.18 




21.17 


21.19 


Voltages (V) 


-29.42 


-29.42 




-29.42 


-29.41 


Error (1 O'") 


3.153 


3.153 




3.153 


3.153 


Operations (10^) 


2.45 


2.97 




1.51 


0.84 



Table 10. Optimization Algorithm Results for Three 10% Length Actuators for Equation 

4.23 Desired Shape 



Run Number 


1* 


2* 


3* 


4 


5 


6 


Initial Actuator 


0.0 


0.3 


0.1 


0.1 


0.0 


0.05 


Locations (m) 


0.1 


0.4 


0.3 


0.3 


0.2 


0.25 




0.2 


0.5 


0.6 


0.5 


0.7 


0.65 




0.3 


0.6 


0.8 


0.7 


0.9 


0.85 


Initial Actuator 


17.07 


151.41 


42.46 


41.55 


19.74 


30.18 


Voltages (V) 


2.89 


-109.40 


4.64 


7.65 


26.31 


17.04 




42.39 


-12.45 


-11.57 


-6.15 


-27.09 


-19.37 




-24.74 


1.26 


-26.61 


-22.63 


3.80 


-21.98 


Optimum 


0.0041 


Did Not 


Did Not 


0.0303 


0.0303 


0.0303 


Actuator 


0.1199 


Converge 


Converge 


0.2273 


0.2273 


0.2273 


Locations (m) 


0.2463 






0.6514 


0.6513 


0.6513 




0.7161 






0.8153 


0.8152 


0.8152 


Optimum 


16.73 






25.76 


25.76 


25.76 


Actuator 


15.69 






20.92 


20.93 


20.92 


Voltages (V) 


13.57 






-16.05 


-16.04 


-16.04 




-29.47 






-21.81 


-21.81 


-21.81 


Error (10‘^) 


27.88 






8.540 


8.540 


8.540 


Operations 

il2!) 


10.7 






6.41 


5.85 


3.38 



Table 11. Optimization Algorithm Results for Four 10% Length Actuators for Equation 

4.23 Desired Shape 



81 



(ui) uieog 



Actual vs. Dccirod Beam Shapes 




Actual 

D«slrcd 



Figure 35. Optimal Locations for Four 10% Length Actuators for Equation 4.23 

Desired Shape 
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Aciiial vs. Dscirdd Beam Shapes 




Aetual 

Desirad 



Figure 36. Actuator Locations Resulting in Local Minimum of Cost Function for 

Equation 4.23 Desired Shape 



F. EFFECT OF NUMBER OF ACTUATORS FOR A GIVEN TOTAL 
ACTUATOR LENGTH 

The optimization algorithm was run for a total length of actuators equal to 50% of 
the beam length divided into one, two, three, four or five actuators to analyze the effect of 
dividing a given total actuator length into increasing numbers of smaller actuators. The 
algorithm was tasked with determining the optimum location for the actuators on a 0.7875 
mm (0.031 in.) thick 7075-T6 aluminum cantilever beam to best approximate a desired 
parabolic beam shape with a tip displacement equal to 0.1% of the beam length. The 
beam curvature per unit actuator input voltage was calculated for a symmetrically 
mounted pair of 0.26 mm (10 mil) thick Navy Type II piezoceramic actuators covering 
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75% of the width of the beam. The results for each number of actuators are tabulated in 
Table 12 and plotted in Figure 37 through Figure 41. These results demonstrate that 
increasing the number of actuators making up the total actuator length contributes 
substantially to minimizing the error in approximating the desired shape. Doubling the 
number of actuators from one to two and from two to four each resulted in about an order 
of magnitude decrease in the error cost function. 



Number of 
Actuators 


1 


2 


3 


4 


5 


Actuator 


0.5 


0.25 


0.1667 


0.125 


0.1 


Length(s) (m) 












Optimum 


0.0852 


0.0600 


0.0450 


0.0359 


0.0298 


Actuator 




0.4684 


0.3426 


0.2646 


0.2155 


Location(s) (m) 






0.6432 


0.4955 


0.4031 










0.7265 


0.5906 












0.7782 


Optimum 


3.429 


3.824 


3.993 


4.087 


4.147 


Actuator 




4.159 


4.355 


4.461 


4.528 


Voltage(s) (V) 






4.350 


4.462 


4.529 










4.456 


4.529 












4.522 


Error (x 10'*^) 


408.9 


54.38 


14.15 


5.166 


2.310 



Table 12. Optimization Algorithm Results for 50% Total Actuator Length for Parabolic 

Desired Shape 



84 



Boam Dkplaeamart Y (m) 




X(m) 

Figure 37. Optimal Location for One 50% Length Actuator 
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B«am Dfeplaeamant Y (m) 



1 

Q.9 
Q.e 
07 
.6 
.5 
0.4 
0.3 
0.2 
0.1 
0 

0 0.2 0.4 0.0 0.e 1 

X(m) 

Figure 38. Optimal Location for Two 25% Length Actuators 



Actual 

Doslrad 



Actual vs. Dcclfdd Baain Shapas 
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B«am Dlaplae«mank Y (m) 



Aetuaivs. Daclrad Beam Shapes 




Actual 

Desired 



X(m) 

Figure 39. Optimal Location for Three 16.67% Length Actuators 
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Beam DkplaeemantY (m] 



Actual vs. Declrdd Beam Shapes 




Actual 

Ddslrad 



X(im) 

Figure 40. Optimal Location for Four 12.5% Length Actuators 
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Beam Dfeplacamont Y (m) 



Actual vs. Ddctrdd Boann Shapes 




Actual 

Ddslrdd 



X(m) 

Figure 41. Optimal Locations for Five 10% Lenph Actuators 
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V. 



EXPERIMENTAL ANALYSIS 



Experimental analysis was conducted in two parts. The first part consisted of 
comparison of the performance of several actuator configurations in producing beam 
curvature with the predictions fi-om the corresponding Euler-Bemoulli models. The 
Euler-Bemoulli predictions neglected the effects of transverse stresses produced by the 
piezoceramic, assumed a linear strain distribution across the thickness of the composite 
actuator/beam and assumed a linear relationship between beam curvature and actuator 
input voltage. The second set of experiments examined the deformation of a beam due to 
four single actuators in the optimal locations to best approximate a parabolic displacement 
profile. 

A. COMPARISON OF ACTUATOR CONFIGURATIONS 

1 . Physical Configuration 

A schematic of the equipment setup for the experiment to compare the 
effectiveness of actuator configurations is shown in Figure 42. Piezoceramic actuators 
were bonded to a 0.031 in. (0.79 mm) thick 7075-T6 aluminum beam which was 
cantilevered such that its length was horizontal and width was vertical to allow bending of 
the beam to take place in the horizontal plane. The single actuator, symmetric pair and 
stacked pair configurations shown in Figure 5, Figure 6 and Figure 7, respectively, were 
each tested on separate beams. Each beam had an overall length of 24 in., 2 in. of which 
was held in the clamp, and a width of 0.875 in. All actuators were 2.5 in. long, 0.75 in. 
wide and 0.26 mm thick Navy Type II piezoceramic. 

A XANALOG nonlinear dynamic systems modeling, simulation and real-time 
controls testing system was used to control actuator input and record experimental data. 
The XANALOG is a PC-based system incorporating an internal digital signal processor 
controlled by NL-SIM software running in a Windows environment. NL-SIM provides a 
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graphical interface which allows a control system model to be constructed using block 
diagrams and then compiled and operated in real time. System inputs and outputs can be 
displayed in both digital and graphical formats while the system is running and recorded to 
data files for post-experiment processing. 




Figure 42. Experimental Setup for Comparison of Actuator Configurations 

Beam tip displacement was measured at a distance of 20. 125 in. from the center of 
each actuator using an NAiS ANL1651AC infrared laser analog displacement sensor with 
an output of 0.1 volts per millimeter of displacement. Beam surface strain was also 
measured directly opposite the center of the single and stacked pair actuator 
configurations in the beam length and width directions using Measurements Group Inc. 
CEA-13-125WT-350 bi-directional strain gages with a 2.16 gage factor. Each gage was 
connected to a Measurements Group Model 2 120 A Strain Gage Conditioner powered by 
a Model 21 10 A Power Supply. The Model 2 120 A output gains were set to 1732 to 
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provide an output of 4 volts per 1000 microstrain measured by each strain gage. Analog 
to digital conversion of the laser and strain gage conditioner outputs was performed by an 
eight-channel CIO-DAS16/330 A/D board and digital to analog conversion of the actuator 
command signal from the XANALOG system was performed by a six-channel CIO- 
DDA06 D/A board. Both boards were manufactured by Computer Boards, Inc. The +10 
volt analog output of the DDA06 was amplified by a 15:1 analog amplifier to provide a 
+150 volt range of input to each actuator. 

The command signal generator shown in Figure 43 was constructed using the NL- 
SIM software. The gain of 1 .66667 was selected to allow integer command inputs to 
produce corresponding multiples of 25 volts input to the piezoceramic actuator when 
amplified by the analog amplifier. The second order filter shown was included to provide 
a gradual response to changes in the command input in order to prevent excitation of 
vibration in the beam. The filter’s damping factor, was set to a value of 1 for critical 
damping and the filter’s natural frequency, co„, was set to 0.1047 rad/s, corresponding to a 
value of 1/60 Hz. The response of the filter allowed the XANALOG output to the analog 
amplifier to reach steady state to four significant figures approximately 1.5 minutes after a 
change in the commanded input was entered. 




Figure 43. Block Diagram of XANALOG Actuator Command Signal Generator 
Implementation for Single Actuator Tests 



2. Experimental Procedure 

All of the experimental equipment was powered for at least 30 minutes prior to 
taking measurements to allow steady state conditions to be reached. The XANALOG 
system was configured to sample commanded actuator input and laser displacement sensor 
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and strain gage output at one second intervals to allow recorded data to be averaged in 
order to minimize the effects of noise in the sensor outputs. The laser and strain gage 
conditioner outputs were adjusted to indicate an average zero value as measured by the 
XANALOG system prior to the start of each measurement run. Each measurement run 
consisted of recording sensor outputs at actuator inputs of 0, 25, 50, 75, 100, 125, 150 
and again 0 volts. Each actuator input was allowed to remain steady for 60 seconds to 
provide 60 sensor output data points for each actuator input to be averaged. The next 
voltage level was commanded immediately after each 60 second measurement and allowed 
two minutes to reach steady state before the next 60 second measurement was initiated. 
At least three measurement runs were conducted for each experiment to verify 
repeatability of the results. Data for each measurement run was saved to a MATLAB 
.MAT data file for processing and plotting. Measurements were taken for the single 
actuator acting in both compression and tension to compare the effectiveness of 
unsymmetric actuators for these two cases. Measurements for the stacked actuator pair 
were taken for the actuator acting in compression, while measurements for the symmetric 
actuator pair were taken with one actuator acting in tension and the other in compression. 

3. Experimental Results 

Measurement data for the single actuator beam is presented in Table 13 and 
plotted in Figure 44 for the actuator acting in compression and Table 14 and Figure 45 for 
the actuator acting in tension. Data for the symmetric actuator and stacked actuator pair 
experiments are presented in Table 15 and Figure 46 and Table 16 and Figure 47, 
respectively. Positive strain indicates tensile strain for the single actuator acting in tension 
and compressive strain for all other experiments. All runs indicated in the tables and 
figures were conducted consecutively with less than five minutes between the completion 
of one run and the start of the next. The first run for all experiments except the single 
actuator acting in compression was the first use of that actuator on the day of the 
experiment. The first measurement run for the single actuator acting in compression was 
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conducted immediately after preliminary testing of the experimental equipment which used 
the same actuator. 



Significant hysteresis, or residual strain after the removal of applied voltage from 
the actuator, was apparent on the first run for all but the single actuator compression 
experiment. Udd [Ref. 14, p. 503] notes that hysteresis occurs in piezoceramic materials 
as the result of “energy dissipation due to internal sliding events in the polycrystaline 
piezoelectric body” and t 3 q)ically ranges from 0.1 to 10% of the strain due to an applied 
electric field prior to the removal of that field. Residual beam tip displacements of 31%, 
23% and 36% of the displacement at maximum input voltage were noted here for the 
single, symmetric and stacked actuator beams respectively. The hysteresis effect is much 
less apparent in the following measurement runs for each actuator since the zero- 
displacement reference was adjusted at the start of each run. Measurements for the single 
actuator acting in compression and runs 2 through 4 of the other actuator experiments are 
therefore displacements and strains above any hysteresis level in the material at the start of 
the run. The displacement measurements for these runs demonstrated excellent 
repeatability for all of the actuators. The repeatability of strain measurements for all but 
the single actuator beam was also excellent. Strain measurements for the single actuator 
showed somewhat more variation from run to run but still reflected the same general trend 
as the other measurements. 

The dashed lines in Figure 44 through Figure 47 depict the linear variation of beam 
tip displacement and beam surface strain predicted by the Euler-Bemoulli model, given by 
Equation 4. 14 for one actuator: 



measurement point at the beam tip. It can be seen from the figures that the Euler- 




(5.1) 



where 




between the center of the actuator and the 



95 



Bernoulli model provides a good approximation of the actual behavior of each 
actuator/beam combination at low actuator input voltage but differs increasingly from the 
actual behavior as voltage is increased. The actual relationship between beam tip 
displacement and actuator input voltage can be seen to be nonlinear. This is consistent 
with the fact that the published piezoelectric coefficient, <i?/, of the piezoceramic 
represents a linear approximation of the material behavior over a specified operating 
range. The dotted lines in Figure 44 through Figure 47 represent second-order least 
squares curve fits of the form: 

y(F) = c„F'+c,F + c, (5.2) 



to the observed displacements as a function of voltage. Coefficients determined for these 
curve fits and predicted values of C/ from Equation 5.1 are given in Table 17. It should be 
noted that the values obtained for are negligibly small in magnitude. An expression for 
a second-order curve fit of actuator curvature as a function of actuator input voltage can 
be obtained by dividing Equation 5.2 by the actuator length // and the measurement 



distance 



A , 

X - X, - — I and omitting the negligible constant term: 
. 2 > 



c.v^+cy 

= -r -j\ KV 






(5.3) 



Values obtained for the ko and ki coefficients and the predicted curvature per unit voltage 
are given in Table 1 8. 

The percentage of actual longitudinal strain in relation to the predicted strain was 
roughly 80% of the percentage of actual beam tip displacement in relation to the predicted 
displacement for the uns 3 anmetric actuator configurations and roughly 140% for the 
symmetric actuator configuration. These percentages are too great to be explained by 
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measurement inaccuracies and variation of the beam and actuator material dimensions, 
which were determined to be within 2% of nominal values. Since strain was measured at 
the surface of the aluminum beam for the unsymmetric configurations and at the surface of 
the piezoceramic for the symmetric configuration, these results suggest that the strain 
distribution was not truly linear through the composite beam cross-section as shown in 
Figure 4. Rather, the strain at the aluminum beam surface was less than predicted by the 
assumed linear strain distribution but greater than predicted at the piezoceramic surface. 
The beam tip displacement measurements nonetheless provide a measure of the overall 
effect of the strain distribution over the entire actuator section of the beam at producing 
curvature of the beam’s neutral axis. 

Measured transverse strain tended to exceed the measured longitudinal strain by 
between 10% and 20%. The longitudinal strain had been expected to be significantly 
greater than the transverse strain due to the significantly larger moment of inertia of the 
beam cross-section about its longitudinal axis than about its transverse axis. The larger 
longitudinal moment of inertia should result in significantly less longitudinal curvature, and 
hence transverse strain, produced by the actuator. The strain gage conditioners for the 
longitudinal and transverse strain gages were swapped to eliminate the possibility of 
measurement equipment contributing to this result but found not to be a causing factor. 
The magnitude of the transverse strain measurements indicate that significant curvature 
about the beam’s longitudinal axis does in fact occur at least locally near the center of the 
actuator. 
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Run 


Predicted 


1 


2 


3 


Tip Displacement - 
mm (% of predicted) 










25 V 


0.51 


0.57(112) 


0.57 (112) 


0.57 (113) 


50 V 


1.02 


1.20(118) 


1.18(116) 


1.19(116) 


75 V 


1.53 


1.92 (125) 


1.89(124) 


1.90(124) 


100 V 


2.04 


2.65 (130) 


2.65 (130) 


2.65 (130) 


125 V 


2.58 


3.47 (134) 


3.47 (134) 


3.46(134) 


150 V 


3.06 


4.38 (143) 


4.37 (143) 


4.32 (141) 


0 V 


0.00 


0.14 


0.08 


0.07 


Longitudinal Strain - 
p£ (% of predicted) 










25 V 


7.8 


6.7 (85) 


9.4 (120) 


8.6 (110) 


50 V 


15.5 


13.0(84) 


17.9 (115) 


15.1 (98) 


75 V 


23.3 


22.2 (95) 


28.2 (121) 


24.9 (107) 


100 V 


31.1 


33.3 (107) 


37.5 ( 120) 


34.6 (111) 


125 V 


38.8 


43.5 (112) 


47.5 (122) 


45.9 (118) 


150 V 


46.6 


56.6(121) 


59.3 (127) 


57.3 (123) 


0 V 


0.0 


2.0 


5.4 


0.3 


Transverse Strain -pc 










25 V 




8.6 


10.3 


9.3 


50 V 




17.9 


20.4 


17.7 


75 V 




28.7 


32.3 


29.3 


100 V 




41.7 


43.5 


41.1 


125 V 




54.0 


56.0 


54.2 


150 V 




69.3 


71.1 


68.9 


0 V 




1.4 


4.9 


0.0 



Table 13. Beam Tip Displacement and Strain Measurements for Single Actuator Acting 

in Compression 
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Figure 44. Beam Tip Displacement and Longitudinal Strain Measurements Single 

Actuator Acting in Compression 
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Run 


Predicted 


1 


2 


3 


4 


Tip Displacement - 
mm (% of predicted) 












25 V 


0.51 


0.87 (171) 


0.56(110) 


0.58(115) 


0.58 (113) 


50 V 


1.02 


1.87 (183) 


1.23 (121) 


1.23 (120) 


1.22 (119) 


75 y 


1.53 


2.92 (191) 


1.95 (127) 


1.94(127) 


1.94 (127) 


100 V 


2.04 


4.03 (198) 


2.72 (134) 


2.71 (133) 


2.70(133) 


125 V 


2.58 


5.18(201) 


3.55 (137) 


3.53 (137) 


3.54 (137) 


150 V 


3.06 


6.33 (207) 


4.50 (147) 


4.45 (145) 


4.42 (144) 


0 V 


0.00 


1.94 


0.14 


0.09 


0.09 


Longitudinal Strain - 
ps (% of predicted) 












25 V 


7.8 


5.2 (67) 


6.6 (85) 


7.1 (91) 


7.8 (100) 


50 V 


15.5 


11.4 (74) 


13.8(89) 


17.2 (111) 


16.7(108) 


75 V 


23.3 


19.3 (83) 


23.2 (100) 


26.0 (112) 


26.2 (113) 


100 V 


31.1 


30.5 (98) 


33.1 (106) 


35.8 (115) 


35.8 (115) 


125 V 


38.8 


45.2 (116) 


42.9 (110) 


46.3 (119) 


46.8 (120) 


150 V 


46.6 


59.2 (127) 


55.3 (119) 


58.2 (125) 


57.5 (123) 


OV 


0.0 


1.0 


-1.1 


1.3 


2.0 


Transverse Strain -pe 












25 V 




14.8 


8.9 


9.2 


9.1 


50 V 




29.7 


18.3 


20.2 


18.8 


75 V 




42.3 


29.8 


30.6 


30.1 


100 V 




63.5 


41.6 


43.0 


41.9 


125 V 




80.9 


54.3 


54.7 


54.4 


150 V 




97.9 


68.6 


69.0 


66.9 


0 V 




30.7 


2.2 


2.5 


2.0 



Table 14. Beam Tip Displacement and Strain Measurements for Single Actuator Acting 

in Tension 
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Figure 45. Beam Tip Displacement and Longitudinal Strain Measurements for Single 

Actuator Acting in Tension 
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Run 


Predicted 


1 


2 


3 


4 


Tip Displacement - 
mm (% of predicted) 












25 V 


0.71 


0.87 (122) 


0.83 (117) 


0.81 (114) 


0.82 (115) 


50 V 


1.43 


1.96(137) 


1.81 (127) 


1.81 (126) 


1.80 (126) 


75 V 


2.14 


3.31 (155) 


2.97 (139) 


2.98 (139) 


2.98 (139) 


100 V 


2.86 


5.00 (175) 


4.20 (147) 


4.18(146) 


4.13 (144) 


125 V 


3.57 


6.71 (188) 


5.37 (150) 


5.28 (148) 


5.23 (146) 


150 V 


4.28 


8.46 (198) 


6.80(158) 


6.67 (156) 


6.57 (153) 


OV 


0.00 


1.93 


0.37 


0.28 


0.21 


Longitudinal Strain - 
|i£ (% of predicted) 












25 V 


14.4 


26.2(182) 


23.4 (163) 


24.7 (172) 


22.3 (155) 


50 V 


28.8 


59.3 (206) 


51.6 (179) 


51.0 (177) 


49.3 (171) 


75 V 


43.1 


100.0 (232) 


80.1 (186) 


81.1 (188) 


79.4 (184) 


100 V 


57.5 


144.9 (252) 


114.8 (200) 


113.5 (197) 


111.7 (194) 


125 V 


71.9 


196.4 (273) 


151.6(211) 


149.6 (208) 


147.1 (205) 


150 V 


86.3 


252.2 (292) 


196.9 (228) 


192.7 (223) 


188.2(218) 


0 V 


0.0 


68.6 


13.2 


10.4 


7.5 


Transverse Strain -p£ 












25 V 




28.0 


24.7 


26.1 


25.0 


50 V 




61.6 


54.9 


54.6 


53.3 


75 V 




102.6 


85.9 


87.1 


85.6 


100 V 




148.9 


122.8 


122.2 


120.7 


125 V 




201.7 


162.3 


160.8 


159.2 


150V 




260.6 


211.1 


207.0 


202.8 


OV 




59.4 


11.6 


10.3 


7.5 



Table 15. Beam Tip Displacement and Strain Measurements for Symmetric Actuator Pair 
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Figure 46. Beam Tip Displacement and Longitudinal Strain Measurements for 

Symmetric Actuator Pair 
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Run 


Predicted 


1 


2 


3 


4 


Tip Displacement - 
mm (% of predicted) 












25 V 


0.56 


1.03 (184) 


0.57 (102) 


0.58(104) 


0.57 (102) 


50 V 


1.12 


2.23 (199) 


1.25 (112) 


1.26(113) 


1.26(112) 


75 V 


1.69 


3.43 (203) 


2.05 (122) 


2.02 (120) 


2.01 (119) 


100 V 


2.25 


4.74(211) 


2.95 (131) 


2.90 (129) 


2.88(128) 


125 V 


2.81 


6.17 (220) 


3.96 (141) 


3.87 (138) 


3.81 (136) 


150 V 


3.37 


7.70 (228) 


5.18(154) 


4.99 (147) 


4.94(147) 


OV 




2.79 


0.35 


0.19 


0.16 


Longitudinal Strain - 
ps (% of predicted) 












25 V 


8.5 


14.1 (166) 


8.2 (97) 


6.2 (73) 


7.2 (70) 


50 V 


17.0 


29.0 (171) 


16.3 (96) 


14.4 (85) 


17.5 (85) 


75 V 


24.5 


45.0(184) 


26.2 (107) 


24.7 (101) 


28.0 (95) 


100 V 


33.9 


65.8 (194) 


37.0 (109) 


34.6 (102) 


40.4 (100) 


125 V 


42.4 


83.2 (196) 


49.2 (116) 


45.9 (108) 


54.7 (108) 


150 V 


50.9 


102.7(202) 


63.9 (126) 


59.8 (118) 


70.4 (115) 


OV 




43.4 


5.7 


2.8 


1.5 


Transverse Strain -pc 












25 V 




18.1 


9.4 


7.8 


7.2 


50 V 




36.9 


19.2 


17.5 


17.5 


75 V 




56.4 


30.6 


28.8 


28.0 


100 V 




76.3 


43.5 


41.2 


40.4 


125 V 




97.8 


58.5 


55.2 


54.7 


150 V 




122.0 


76.4 


71.8 


70.4 


OV 




51.2 


6.4 


2.9 


1.5 



Table 16. Beam Tip Displacement and Strain Measurements for Stacked Actuator Pair 

Acting in Compression 
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Figure 47. Beam Tip Displacement and Longitudinal Strain Measurements for Stacked 

Actuator Pair Acting in Compression 
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Beam 


Predicted Ci 
(10‘VnW) 


Co 

(10"’ mm/V^) 


Cl 

(10'^ mmA'^) 


C2 

(mm) 


Single - Compression 


2.04 


5.05 


2.15 


0.00 


Single - Tension 


2.04 


5.15 


2.20 


0.00 


Symmetric Pair 


2.86 


6.79 


3.46 


-0.04 


Stacked Pair 


2.25 


8.72 


2.03 


0.01 



Table 17. Coefficients of Curve Fits for Observed Beam Tip Displacement as a Function 

of Actuator Input Voltage 



Beam 


Predicted ki = Ki 


ko 


ki 




(lO"* rad/mA^) 


(10-^ rad/mA^^) 


(10'^ rad/mA^) 


Single - Compression 


6.28 


1.56 


6.61 


Single - Tension 


6.28 


1.58 


6.77 


Symmetric Pair 


8.80 


2.09 


10.66 


Stacked Pair 


6.93 


2.69 


6.26 



Table 18. Coefficients of Curve Fits for Observed AcWator Curvature as a Function of 

Actuator Input Voltage 



B. BEAM SHAPE CONTROL WITH FOUR ACTUATORS 



1. Physical Configuration 

A schematic of the equipment setup for the four-actuator shape control experiment 
to compare the effectiveness of actuator configurations is shown in Figure 48. Eight 
piezoceramic actuator patches were bonded in four groups of two to one side of a 0.03 1 
in. (0.79 mm) thick 7075-T6 aluminum beam which was cantilevered such that its length 
was horizontal and width was vertical to allow bending of the beam to take place in the 
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horizontal plane. The beam was also supported at approximately two-thirds of its length 
by an air pad riding on a granite table. The beam had an overall length of 45 in., 2 in. of 
which was held in the clamp, and a width of 1.625 in. All actuator patches were 2.5 in. 
long, 1.50 in wide and 0.26 mm thick Navy Type II piezoceramic, hence each of the four 
actuator groups had a total length of 5 in. The actuators were placed at locations on the 
beam shown in Table 19 which were determined by the optimization algorithm to best 
approximate a parabolic deflected beam shape with four 5 in. long actuators on a 43 in. 
long beam. The XANALOG system was again used to control actuator input and record 
experimental data and the ANL1651AC laser displacement sensor was used to measure 
displacement at selected points along the length of the beam. Beam surface strain was 
also measured directly opposite the center of the actuator patch nearer to the cantilevered 
end of the beam in each of the four actuator groups using a Measurements Group Inc. 
CEA-13-125WT-350 bi-directional strain gage opposite the first actuator and a CEA-13- 
125UW-350 unidirectional strain gage opposite the other three actuators. Strain gage 
output was obtained with a Gould Model 56-1301-00 DC/Bridge/Transducer Signal 
Conditioner for each gage. DAS16/330 and DDA06 boards were used for A/D and D/A 
conversion, respectively. Each ±10 volt analog output channel of the DDA06 was 
amplified by a 15; 1 analog amplifier to provide a +150 volt range of input to each 
actuator. 



Actuator 


Distance from 


Percentage of Actuator 3 




Cantilevered End 


Curvature for Parabolic 




(in/cm) 


Shape Approximation 


1 


1.70/4.31 


91.4 


2 


11.52/29.26 


100.0 


3 


21.43/54.43 


100.0 


4 


31.34/79.60 


99.9 



Table 19. Actuator Locations and Curvature Ratios for Parabolic Shape Approximation 
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Figure 48. Setup of Four Actuator Shape Control Experiment 



1. Experimental Procedure 

A series of measurements was initially taken to characterize the relationship 
between actuator curvature and input voltage for the single actuator configuration on the 
beam used for this set of experiments. Beam tip displacement measurements were taken 
with voltage applied to the first actuator only and then with voltage apphed to the second 
actuator only using the procedures described for the previous set of experiments. The 
beam tip displacement measirrement point was 38 in. from the center of the first actuator 
and 28.13 in. from the center of the second actuator. 

Two shape control experiments were then performed. Each experiment consisted 
of applying voltages predicted to produce the optimum approximation of a parabolic 
deformation profile along the beam’s length from its initial undeformed shape. Voltages 
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for each experiment were calculated to produce curvature of the four actuators in the 
proportions listed in Table 19. Identical voltages were applied to actuators 2 through 4 
since there was negligible difference between the curvatures required from these actuators. 
Corresponding input voltages for actuator 1 were determined for the first experiment by 
assuming a linear relationship between curvature and voltage as developed in Chapter II. 
Actuator I voltages for the second experiment were determined by solving Equation 5.3 
for the curvature produced by each actuator 2-4 input voltage and then solving Equation 
5.3 for the voltage to produce actuator 1 curvature equal to 91.4% of the actuator 2-4 
curvature. Coefficients ko and ki for Equation 5.3 of 5.19x10'^ rad/mA^^ and 6.26x10''’ 
rad/mA'^, respectively, were obtained from a curve fit of beam tip displacement to actuator 
1 input voltage as described in the following Experimental Results section. 

Beam displacement was measured using the laser displacement sensor at locations 
of 18 in., 30 in. and 42.22 in. from the cantilevered end of the beam for each of the two 
experiments. It should be noted that the laser measurement axis was oriented 
perpendicular to the longitudinal axis of the undeformed beam. Increasing beam 
displacement therefore resulted in the displacement measurement points slightly further 
down the length of the beam than the nominal measurement locations. However, the small 
magnitude of the measured displacements relative to the length of the beam for these 
experiments made this effect insignificant. 

Three sets of measurements were taken for each of the two experiments at each 
measurement point. Each measurement run consisted of recording sensor outputs at 
actuator inputs given for the input stages in Table 20 and Table 21 for the first and second 
experiment, respectively. Each actuator input was allowed to remain steady for 60 
seconds to provide 60 sensor output data points for each actuator input to be averaged. 
The next input voltage stage was commanded immediately after each 60 second 
measurement and allowed two minutes to reach steady state before the next 60 second 
measurement was initiated. Data for each measurement run was saved to a MATLAB 
.MAT data file for processing and plotting. All measurements were taken with the 
actuators acting in compression. The actuators were cycled through the input voltages 
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corresponding to one measurement run prior to the recording of data to minimize the 
effects of hysteresis on the repeatability of the measurements. 



Input Stage 


Actuator 1 
Input Voltage (V) 


Actuator 2,3,4 
Input Voltage (V) 


1 


0 


0 


2 


45.7 


50 


3 


91.4 


100 


4 


137.2 


150 


5 


0 


0 



Table 20. Linear Model Actuator Input Voltages for Shape Control Experiment I 



Input Stage 


Actuator 1 
Input Voltage (V) 


Actuator 2,3,4 
Input Voltage (V) 


1 


0 


0 


2 


46.5 


50 


3 


93.6 


100 


4 


141.1 


150 


5 


0 


0 



Table 21. Nonlinear Model Actuator Input Voltages for Shape Control Experiment 2 



3. Experimental Results 

Displacement and strain measurement data for beam deformation due to the first 
actuator is presented in Table 22 and plotted in Figure 50. Data for beam deformation due 
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to the second actuator is presented Table 23 and Figure 52. Positive strain values indicate 
compressive strain. All runs indicated in the tables and figures were conducted 
consecutively with less than five minutes between the completion of one run and the start 
of the next. The first run for each actuator was the first use of that actuator on the day of 
the experiment. Some hysteresis was observed on all measurement runs but most apparent 
on the first run for each actuator, with residual beam tip displacements of 16% and 24% of 
the displacement at maximum input voltage occurring on the first run for the first and 
second actuator, respectively. Repeatability of the displacement measurements for both 
actuators and strain measurements for the first actuator was good for runs 2 through 4, 
but the strain measurements for the second actuator were not as reliable. Strain 
measurements for runs 1 and 3 of the second actuator were substantially lower than 
corresponding measurements for other runs of the first and second actuator. The final 
zero input strain measurements for both of these runs indicated residual tensile rather than 
compressive strain, as denoted by the asterisks in Table 23. The combination of low 
measurement values and final tensile strain measurements suggest that either a shift in the 
strain gage conditioner zero strain output voltage occurred shortly after the start of each 
of these runs or the bond between the actuator 2 strain gage and the beam may have been 
poor, leading to inconsistent results. 

The dashed lines in Figure 50 and Figure 52 represent the relationship between 
actuator input voltage and beam tip displacement predicted by the linear Euler-Bemoulli 
model, while the dotted lines represent second order least squares curve fits of the data 
from runs 2 through 4 for each actuator. Coefficients of Equation 5.2 for each of the 
curves are listed in Table 24 and the corresponding Equation 5.3 coefficients obtained are 
listed in Table 25. The ko and ki coefficients obtained for the two actuators agreed to 
within 6% and 2%, respectively. 

It can be seen from the figures that the Euler-Bemoulli model again provides a 
reasonable approximation of the behavior of the beam as a fimction of input voltage, as 
was observed for the previous beams tested. However, unlike the previous experiments, 
the Euler-Bemoulli model overpredicted the tip displacement of the 1.625 in. wide beam 
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for actuator input voltages below about 75 volts and underpredicted the tip displacements 
at higher voltages. Comparison of the percentages of predicted displacement obtained for 
the 1.625 in. wide beam in Table 22 and Table 23 with the corresponding percentages of 
predicted displacement obtained for the 0.875 in. wide beam in Table 13 and Table 14 
reveals that the percentage of predicted displacement obtained was always less for the 
wider beam. The reason for this is that the actuator produces curvature of the beam about 
both its transverse and longitudinal axes, a fact supported by the significant transverse 
strain measurements obtained. Curvature about the longitudinal axis increases the 
moment of inertia of the beam’s cross-section, as shown in Figure 51, resulting in a 
greater resistance of the beam to bending about its transverse axis. The increase in 
moment of inertia is more significant for the wider beam than for the narrower beam. 

The results of shape control experiment 1 are tabulated in Table 26 and plotted in 
Figure 52 and Figure 53. The “linear model predicted” column of Table 26 and the 
dashed lines on the figures represent the desired displacement of the beam along its length 
to approximate a parabolic shape predicted by the Euler-Bemoulli model of actuator 
curvature and Equations 4.10 through 4.14. The “nonlinear model predicted” column in 
the table and dotted lines on the figures represent the displacement of the beam predicted 
by Equations 4.10 through 4.14 with the curvature terms KiVj replaced by Equation 5.3 
with the empirical coefficients ko and ki equal to those obtained from the experiment to 
characterize the first actuator given in Table 25. The measured displacements displayed 
excellent repeatability and can be seen to agree much more closely with the beam shape 
profile predicted using the nonlinear empirical model of actuator curvature than for the 
profile predicted by the linear model. As for the single actuator tests, the linear model 
overestimates the actuator performance for the low input voltages and underestimates 
actuator performance at higher voltages. 

The results of shape control experiment 2 are presented in Table 27, Figure 54 and 
Figure 55. The “nonlinear model predicted” column in the table and dotted lines on the 
figures represent the desired displacement of the beam to approximate a parabolic shape 
determined using Equations 4.10 through 4.14 with the curvature terms replaced by 
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Equation 5.3, with the empirical coefficients ko and ki equal to those for the first actuator 
in Table 25. The “linear model predicted” column in the table and dashed lines on the 
figures represent the displacement of the beam predicted by the linear Euler-Bemoulli 
model of actuator curvature and Equations 4.10 through 4.14 for the specified input 
voltages. The measured displacements at the 18 in. and 30 in. measurement points 
displayed excellent repeatability, but the displacements at the 42.22 in. measurement point 
showed somewhat more variation. The latter measurements were the first series taken for 
this experiment and therefore the most affected by hysteresis. Residual displacements of 
about 6% of the displacement at maximum input voltage were observed for runs 1 and 2 
at the 42.22 in. measurement point, but less than 2% residual displacement was observed 
on run 3 and for the measurements at the other measurement points. With the exception 
of runs 1 and 2 at the 42.22 in. measurement point, the measured displacements displayed 
much better agreement with the desired beam shape profile predicted using the nonlinear 
empirical model of actuator curvature than with the profile predicted by the linear model. 

The results of the shape control experiments demonstrate the importance of 
considering the nonlinear behavior of the piezoceramic in determining the appropriate 
actuator inputs for precision shape control applications. The significant effects of 
hysteresis observed indicate that this effect must also be considered. Actual shape control 
applications should incorporate some form of feedback of the beam’s shape to ensure that 
the desired deformation profile is attained. 
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Run 


Predicted 


1 


2 


3 


4 


Tip Displacement - 
mm (% of predicted) 












25 V 


1.98 


2.27(115) 


1.13 (57) 


1.33 (67) 


1.26(64) 


50 V 


3.96 


4.60(116) 


3.20 (81) 


3.49 (88) 


3.52 (89) 


75 V 


5.94 


7.33 (123) 


5.97 (100) 


5.71 (96) 


5.69 (96) 


100 V 


7.92 


10.12(128) 


8.89(112) 


8.71 (110) 


8.19(103) 


125 V 


9.90 


13.33 (135) 


11.55(117) 


11.81 (119) 


10.84(110) 


150 V 


11.88 


17.54 (148) 


15.21 (128) 


15.59(131) 


14.73 (124) 


OV 


0.00 


2.78 


0.93 


1.29 


0.19 


Longitudinal Strain - 
pe (% of predicted) 












25 V 


8.0 


7.2 (90) 


5.5 (69) 


5.4 (68) 


6.5 (82) 


50 V 


16.0 


14.5 (91) 


13.7 (85) 


12.7 (80) 


13.3 (83) 


75 V 


24.0 


23.4 (97) 


21.2 (88) 


21.2 (89) 


21.5(90) 


100 V 


32.0 


32.1 (100) 


30.3 (95) 


29.4 (92) 


30.4 (95) 


125 V 


40.0 


42.2 (106) 


40.1 (100) 


39.2 (98) 


39.5 (99) 


150 V 


48.0 


54.1 (113) 


51.5 (106) 


50.3 (105) 


49.8 (104) 


OV 


0.0 


4.6 


1.8 


2.4 


0.9 


Transverse Strain -ps 












25 V 




8.2 


8.2 


8.0 


8.1 


50 V 




17.6 


17.1 


17.0 


16.3 


75 V 




27.5 


26.8 


26.8 


26.5 


100 V 




38.8 


37.5 


37.3 


36.6 


125 V 




51.6 


49.7 


49.8 


48.7 


150 V 




67.0 


63.7 


63.4 


60.8 


OV 




5.7 


3.4 


3.2 


-0.6 



Table 22. Beam Tip Displacement and Strain Measurements for First Actuator Acting in 

Compression 
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Figure 49, Beam Tip Displacement and Longitudinal Strain Measurements for First 

Actuator Acting in Compression 
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Run 


Predicted 


1 


2 


3 


4 


Tip Displacement - 
mm (% of predicted) 












25 V 


1.47 


1.29(88) 


1.31 (89) 


1.12(76) 


1.27 (86) 


50 V 


2.93 


3.06 (104) 


2.71 (93) 


2.67 (91) 


2.76 (94) 


75 V 


4.40 


5.45 (124) 


4.54 (103) 


4.38 (100) 


4.45 (101) 


100 V 


5.86 


8.08 (138) 


6.44 (110) 


6.41 (109) 


6.22 (106) 


125 V 


7.33 


11.22 (153) 


8.88 (121) 


8.47 (116) 


8.81 (120) 


150 V 


8.79 


15.05 (171) 


11.81 (134) 


11.40 (130) 


11.62 (132) 


OV 


0.00 


3.57 


0.81 


0.31 


0.85 


Longitudinal Strain - 
|iE (% of predicted) 












25 V 


8.0 


3.3 (41) 


5.7 (72) 


4.9 (62) 


5.9(74) 


50 V 


16.0 


7.9 (49) 


11.8 (74) 


10.3 (64) 


12.2 (76) 


75 V 


24.0 


16.6 (69) 


18.9(79) 


16.5 (69) 


19.6 (82) 


100 V 


32.0 


25.6 (80) 


27.2 (85) 


23.6 (74) 


27.0 (84) 


125 V 


40.0 


34.1 (85) 


36.9 (92) 


31.2 (78) 


36.4 (91) 


150 V 


48.0 


44.7 (93) 


47.3 (99) 


40.2 (84) 


45.6 (95) 


OV 


0.0 


-3.6* 


0.8 


-5.5* 


0.8 



Table 23. Beam Tip Displacement and Strain Measurements for Second Actuator Acting 

in Compression 
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Figure 50. Beam Tip Displacement and Longitudinal Strain Measurements for Second 

Actuator Acting in Compression 
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Actuator 


Predicted Ci 


Co 


Cl 


C2 




(lO'^mm/V) 


(10-' mmA^^) 


(10'^ mm/V) 


(mm) 


First 


7.92 


3.18 


5.38 


-0.10 


Second 


5.86 


2.48 


3.91 


0.06 



Table 24. Coefficients of Curve Fits for Observed Beam Tip Displacement as a Function 

of Actuator Input Voltage 



Actuator 


Predicted ki = Ki 


ko 


ki 




(lO'^' rad/mA^) 


(10-^ rad/mA^^) 


(lO'^rad/mA^) 


First 


6.46 


2.59 


4.39 


Second 


6.46 


2.74 


4.31 



Table 25. Coefficients of Curve Fits for Observed Actuator Curvature as a Function of 

Actuator Input Voltage 




Figure 51. Graphic Depiction of Deformation of Beam Cross-Section due to 
Piezoceramic Actuator Acting in Compression 
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Beam Displacement at 
Measurement Locations - mm 


Linear 

Model 

Predicted 


Nonlinear 

Model 

Predicted 


Run 1 


Run 2 


Run 3 


45.7/50/50/50 Volts Input 












X = 18.0 in 


1.73 


1.50 


1.47 


1.25 


1.32 


X = 30.0 in 


4.76 


4.14 


3.97 


3.91 


4.00 


X = 42.2 in 


9.33 


8.16 


8.45 


8.20 


8.12 


91.4/100/100/100 Volts Input 












X = 18.0 in 


3.46 


3.65 


3.57 


3.40 


3.42 


X = 30.0 in 


9.51 


10.11 


9.86 


9.88 


9.91 


X = 42.2 in 


18.67 


19.94 


20.63 


20.21 


20.36 


137.2/150/150/150 
Volts Input 












X = 18.0 in 


5.19 


6.45 


6.03 


5.88 


5.86 


X = 30.0 in 


14.27 


17.91 


17.20 


17.07 


17.07 


X = 42.2 in 


28.02 


35.35 


35.18 


34.73 


34.82 



Table 26. Beam Displacement Measurements for Shape Control Experiment 1 (Linear 

Model Voltages) 
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Actuator Inputs; 4S7/5Q/SQ/SQ Volts 




— 


Lkidar 


• 


Nonibidar 


X 


Run 1 


+ 


Run 2 


o 
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Actuator Inputs: 91 4/1QQ/1QQ/1QQ Volts 




— 
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• 
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X 
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Figure 52, Beam Displacement Measurements for Shape Control Experiment 1 (Linear 

Model Voltages) (1,2 of 3) 
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Bfiam Dfeplac^mant (mm) 



Aetuaioi Inputs; 137.2/1SQ/1SQ/1SQ Volts 




Figure 53. Beam Displacement Measurements for Shape Control Experiment I (Linear 

Model Voltages) (3 of 3) 
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Beam Displacement at 
Measurement Locations - mm 


Linear 

Model 

Predicted 


Nonlinear 

Model 

Predicted 


Run 1 


Run 2 


Run 3 


46.5/50/50/50 Volts Input 












x= 18.0 in 


1.75 


1.52 


1.37 


1.25 


1.30 


X = 30.0 in 


4.80 


4.19 


4.10 


3.85 


3.87 


X = 42.2 in 


9.40 


8.22 


8.99 


9.46 


8.26 


93.6/100/100/100 Volts Input 












X = 18.0 in 


3.52 


3.74 


3.60 


3.48 


3.49 


X = 30.0 in 


9.63 


10.28 


10.10 


9.85 


9.66 


X = 42.2 in 


18.85 


20.18 


22.18 


22.04 


20.27 


141.1/150/150/150 
Volts Input 












X = 18.0 in 


5.30 


6.65 


6.23 


6.16 


6.17 


X = 30.0 in 


14.48 


18.29 


17.76 


17.54 


17.08 


X = 42.2 in 


28.33 


35.91 


37.98 


37.54 


35.45 



Table 27. Beam Displacement Measurements for Shape Control Experiement 2 

(Nonlinear Model Voltages) 
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Beam Dlspbcamant (mm) Baam Dteplacamant (mm) 



Actuator Inputs: 4fi.S.€Q/SQ/SQ Volts 




— 


Ltodai 


• 


Nonitodar 


X 


Run 1 


+ 


Run 2 


o 


Run 3 



Actuator inputs; 03.fi/1QQ/1QQ/1QQ Volts 




Ltodar 

Nonikidar 
X Run 1 
+ Run 2 
o Run 3 



Figure 54, Beam Displacement Measurements for Shape Control Experiment 2 
(Nonlinear Model Voltages) (1,2 of 2) 
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Aeiu&toi Inputs; UI.I/ISQ/tSQ/tSQ Volts 




— 
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Figure 55. Beam Displacement Measurements for Shape Control Experiment 2 
(Nonlinear Model Voltages) (3 of 3) 
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VI. CONCLUSIONS 



Piezoceramic actuators have been shown to provide an effective means of 
controlling the shape of a thin flexible beamlike structure. Experimental analysis has 
demonstrated that the Euler-Bemoulli model provides a good measure of the relative 
effectiveness of various actuator configurations and a reasonable model of actuator 
performance at low input voltage. However, the true nonlinear behavior of piezoceramic 
actuators and the effects of hysteresis and transverse stresses not accounted for by the 
Euler-Bemoulli model have been shown to have significant ramifications for the use of 
piezoceramics in shape control. Hysteresis resulted in residual beam deformation in excess 
of 20% of full deformation after the removal of the initial actuator input voltage for the 
Navy Type II piezoceramics tested, but following displacement measurements were found 
to be highly repeatable. Nonlinear voltage-strain behavior of the piezoceramic produced 
significantly higher displacements than predicted by the Euler-Bemoulli model at higher 
actuator input voltages. Transverse stresses were found to result in reduced effectiveness 
for a given thickness of actuator on wider beams due to the resulting increase in the 
moment of inertia of the beam cross-section. All of these effects result in variation of the 
relationship between actuator curvature and input voltage from the linear behavior 
predicted by the Euler-Bemoulli model. Precision shape control applications must 
therefore consider these effects and should employ a control system incorporating some 
form of strain and/or position feedback to ensure the desired deformation of the beam is 
achieved. 

The amount of beam deformation available from a piezoceramic actuator has been 
shown to depend on the dimensions and material properties of the structure and 
piezoceramic. High piezoceramic elastic modulus, piezoelectric coefficient and operating 
voltage are all desired actuator properties. The optimum total actuator thickness for a 
given beam/actuator configuration has been shown to be equal to a fixed percentage of the 
beam thickness, while the total beam curvature attainable for the optimum beam/actuator 
thickness ratio has been shown to be inversely related to the thickness of the beam. 
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The location of actuators on a beam structure is crucial to determining how closely 
a desired deformation profile can be approximated. Simultaneous optimization of the 
locations and input voltages of a fixed set of actuators to achieve a desired beam 
deformation profile proved to be unreliable due to the differences in the order of the 
actuator location and voltage terms in the optimization cost function. Using embedded 
Nelder and Mead simplex algorithms to separately optimize actuator locations and input 
voltages was found to produce much more reliable results, converging to the same 
optimum solution for a variety of initial conditions. Optimization of locations for different 
numbers of actuators making up a fixed total actuator length has shown that dividing a 
given actuator length into multiple small actuators allows a desired beam shape to be 
approximated more closely than can be achieved with a lesser number of longer actuators. 

Certain sets of initial conditions were found to cause the optimization algorithm to 
fail to converge or converge to an erroneous result. The algorithm was also found to be 
susceptible to converging to a local minimum of the cost function, when one existed and 
initial conditions near the local minimum were selected. Consideration must therefore be 
given to the presence of inflection points in the desired shape function, which can result in 
local minima of the cost function, when selecting initial actuator locations for the 
optimization algorithm. Knowledge of the approximate values of the optimum actuator 
locations should be used in selecting initial conditions whenever possible. Multiple runs of 
the optimization program should be performed using different sets of initial conditions to 
ensure that the globally optimum actuator locations are attained. 

Recommendations for further research include improving the mathematical model 
of the beam to consider transverse effects, incorporating constraints and techniques to 
avoid local minima in the optimization algorithm and exploring the use of feedback to 
control the shape of the beam. A plate theory model can be used to consider the effects of 
the curvature about the beam’s longitudinal axis caused by the actuator acting in both the 
longitudinal and transverse directions for beam’s with high width-to-thickness ratios. 
Modifying the optimization algorithm to include constraints would allow for the possibility 
of desired shape functions which can be most closely approximated by placing actuators 
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end to end or by actuators providing the maximum possible curvature. Techniques such 
as genetic algorithms or simulated annealing could also be included in the optimization 
algorithm to improve the probability of attaining the global minimum rather than a local 
minimum of the cost function. Multiple strain and/or displacement sensors could be 
utilized to provide feedback of the beam’s shape to control inputs to the piezoceramic 
actuators. Sensors such as the laser displacement sensor used in the experiments 
conducted for this thesis and state-of-the-art fiber optic strain sensors are strong 
candidates for this application. 
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APPENDIX 



% PROGRAM TO COMPARE PERFORMANCE OF PIEZOCERAMIC ACTUATOR 
% CONFIGURATIONS 
% 

% This program calculates and plots the following data: 

% 1) Beam curvature per unit voltage 
% 2) Beam curvature at 20 Volts/mil of piezo thickness 
% 3) Beam surface strain per unit voltage (Side of beam 
% opposite actuator for unsymmetric configurations) 

% 4) Neutral axis distance from center of beam 
% versus actuator thickness or actuator width for piezoceramic 
% actuators bonded to a cantilever beam, given input properties 
% of the beam, actuators and the bond between the actuator and 
% the beam. 

% 

% The following actuator configurations are evaluated; 

% 1 ) A single actuator bonded to one side of the beam 
% 2) A single actuator bonded to each side of the beam 
% 3) A dual (stacked) actuator bonded to one side of the beam 
% 4) A dual (stacked) actuator bonded to each side of the beam 

% Set number of points for plots 
Numpts = 100; 

% Set properties of beam 
% Material name 
Beammatl = '7075-T6 Aluminum'; 

% Width (m) 

w = 0.0381; % 1.5 inches 

% Thickness 
t = 0.001575; 

% Elastic modulus (N/m''2) 

E = 7.2el0; 

% Set properties of PZT actuator 
% Material name 
Piezonam = ‘DOD Type II PZT'; 

% Constant width (m) for plots vs. actuator thickness 
wp = 0.0381; % 1.5 inches 

% Width range (m) for plots vs. actuator width 
wprange = linspace(0. 1 , 1 .0,Numpts)* w; 

% Constant thickness (m) for plots vs. actuator width 
tp = 0.00026; 



129 



% Thickness range (m) for plots vs. actuator thickness 
tprange = linspace(0.0001,0.002,Numpts); 

% Elastic modulus (N/m''2) 

Ep = 6.1el0; 

% Piezoelectric coefficient (nv^) 
d31 = 1.71e-10; % DOD Type II PZT 

% Set properties of bond 
% Width (m) 
wb = wp; 

% Thickness (m) 
tb = 0; 

% Elastic modulus (N/m^2) 

Eb = 1.78e9; 

% Set electric field limitation of piezoceramic (volts/mil) 

Volt_mil = 20; % Volts/mil 

% Calculate electric field limitation of piezoceramic (volts/m) 

Vjlimit = Volt_mil/0.0000254; % Volts/m = Volts/miy(0.0000254m/mil) 

% Allow plots to be generated vs actuator width or thickness 

Plottype = input('Generate plots vs actuator thickness (t) or width (w)? Vs'); 

% Allow actuator thickness to be plotted as % of beam thickness 
% or actual thickness for plots vs actuator thickness 
if Plottype ~= 'w' 

Axistype = input('Plot actuator thickness as % of beam thickness (p) or actual thickness 
(a)? Vs'); 
end 

% Create actuator width or thickness range in accordance 
% with plot type 
if Plottype = 'w' 
w = w*ones(l,Numpts); 
wp = wprange; 
wb = wp; 
else 

t = t*ones(l,Numpts); 
tb = tb*ones(l,Numpts); 
tp = tprange; 
end 

% Calculate parameters for Crawley/ Anderson formulation 
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% Single Actuator 

EAl = E*w.*t + Eb*wb.*tb + Ep*wp.*tp; 

ESI =Eb*wb.*tb.*(t/2+tb/2) + Ep*wp.*tp.*(t/2+tb+tp/2); 

Ell = E*w.*t.^3/12 + Eb*wb.*tb.*((t/2).^2 + (t/2).*tb + tb.^2/3) + ... 

Ep*wp*tp.*((t/2+tb).^2 + (t/2+tb).*tp + tp.'^2/3); 

PI = Ep*d31*wp; 

Ml = Ep*d3 1 *wp.*(t/2 + tb + tp/2); 

EpsO_l = (P1.*EI1-M1.*ES1)./(EA1.*EI1-ES1.^2); 

K_V1 =M1./EI1 -ES1./EI1.*(P1.*EI1 - M1.*ES1)./(EA1.*EI1 - ESl.^2); 
yO_l =EpsO_l./K_Vl; 

% Sjonmetric Pair of Actuators 

EA2 = E*w.*t + 2*Eb*wb.*tb + 2*Ep*wp.*tp; 

EI2 = E*w.*t.^3/12 + 2*Eb*wb.*tb.*((t/2).^2 + (t/2).*tb + tb.^2/3) + ... 

2*Ep*wp.*tp.*((t/2+tb).^2 + (t/2+tb).*tp + tp.''2/3); 

M2 = 2*Ep*d3 1 *wp.*(t/2 + tb + tp/2); 

K_V2 = M2./EI2; 

% Stacked Pair of Actuators 

EAl 2 = E*w.*t + 2*Eb*wb.*tb + 2*Ep*wp.*tp; 

ES12 = Eb*wb.*tb.*(t+2*tb + tp) + Ep*wp.*tp.*(t+3*tb+2*tp); 

E112 - E*w.*t.^3/12 + Eb*wb.*tb.*((t/2).^2 + (t/2+tb+tp).^2 + ... 
(t+5/3*tb+tp).*tb) + Ep*wp.*tp.*((t/2+tb).^2 + ... 

(t/2+2*tb+tp).^2 + (t+3*tb+5/3*tp).*tp); 

P12 = 2*Ep*d31*wp; 

M12 = Ep*d3 1 *wp.*(t/2 + tb + tp/2 + t/2 + 2*tb + 3/2*tp); 

Eps0_12 = (P12.*E112-M12.*ES12)./(EA12.*E112-ES12.^2); 

K_V12 = M12./E112 - ES12./E112.*(P12.*E112 - M12.*ES12)./ ... 

(EA12.*E112-ES 12.^2); 
y0_12 =EpsO_12./K_V12; 

% Symmetric Pair of Stacked Actuators 
EA22 = E*w.*t + 4*Eb*wb.*tb + 4*Ep*wp.*tp; 

E122 = E*w.*t.^3/12 + 2*Eb*wb.*tb.*((t/2).^2 + (t/2+tb+tp).^2 + ... 
(t+5/3*tb+tp).*tb) + 2*Ep*wp.*tp.*((t/2+tb).^2 + ... 

(t/2+2*tb+tp).'^2 + (t+3*tb+5/3*tp).*tp); 

M22 = 2*Ep*d3 1 *wp.*(t/2 + tb + tp/2 + 1/2 + 2*tb + 3/2*tp); 

K_V22 = M22./E122; 

% Calculate strain per unit actuator input voltage at beam surface 
% for each actuator configuration 
StrainVl = EpsO_l+K_Vl.*t/2; % Single Actuator 
StrainV2 = K_V2. *t/2; % Symmetric Pair of Actuators 

StrainV 12 = EpsO_ 1 2+K_V 1 2 . *t/2 ; % Stacked Pair of Actuators 
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Strain V22 = K_V22.*t/2; 



% Symmetric Pair of Stacked Actuators 



% Plot results 
figure 

if Plottype == V 

plot(wp./w* 100,K_Vi;b-',wp./w* 100,K_V12;m-.>p./w* 100,K_V2;r-', ... 
wp./w*100,K_V22,’g.') 

text(0.15,0.6,[num2str(tp*1000) ' mm Thick Actuator/' num2str(tb(l)*1000) .. 

' mm Bond Thickness'], 'Units','Normalized') 
xlabelC Actuator Width (%)') 
elseif Axistype — 'p' 

plot(tp./t* 1 00,K_V 1 ,'b-',tp./t* 1 00,K_V1 2,'m-.',tp./t* 1 00,K_V2,'r-', . . . 
tp./t*100,K_V22,'g.') 

text(0.15,0.6,[num2str(wp/w*100) '% Width Actuator/' num2str(tb(l)*1000) . 

' mm Bond Thickness'],Units','Normalized') 
xlabelC Actuator Thickness (% of Beam Thickness)') 
else 

plot(tp* 1 000, K_V 1 ,'b-',tp* 1 000,K_V 1 2,'m-.',tp* 1 000,K_V2,'r-', . . . 
tp*1000,K_V22,'g.') 

text(0.15,0.6,[num2str(wp/w*100) '% Width Actuator/' num2str(tb(l)*1000) . 

' mm Bond Thickness'],'Units','Normalized') 
xlabel('Actuator Thickness (mm)') 
end 

ylabel('(rad/m/V)') 

titleCBeam Curvature per Unit Voltage') 
legend('Single','l Stacked','Sym. Patr','Sym. Stacked',- 1) 
text(0. 1 5,0.9,Piezonam,'Units','Normalized') 

text(0.15,0.8,['Ep = ' num2str(Ep) ' N/m''^2 d31 = ' num2str(d31) ' m/V'], ... 
'Units','Normalized') 

text(0.15,0.7,[num2str(t(l)*1000) ' mm Thick ' Beammatl ' Beam'],'Units', ... 

'Normalized') 

grid 

figure 

if Plottype = 'w' 

plot(wp./w* 1 00, V_tlimit*tp. *K_V 1 ,'b— ',wp./w* 1 00, V_tlimit*tp. *K_V 1 2,'m- . 

wp./w* 100,V_tlimit*tp.*K_V2,'r-',wp./w* 1 00,V_tlimit*tp.*K_V22,'g.') 
text(0.15,0.6,[num2str(tp*1000) ' mm Thick Actuator/" num2str(tb(l)*1000) . 

' mm Bond Thickness'],'Units','Normalized') 
xlabel('Actuator Width (%)') 
elseif Axistype = 'p' 

plot(tp./t* 1 00,V_tlimit*tp. *K_V1 ,'b-',tp./t* 1 00,V_tlimit*tp. *K_V1 2,'m-.', ... 

tp./t* 1 00,V_tlimit*tp. *K_V2,'r-',tp./t* 1 00,V_tlimit*tp. *K_V22,'g.') 
text(0.15,0.6,[num2str(wp/w*100) '% Width Actuator/' num2str(tb(l)*1000) 
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' mm Bond Thickness'],XJnits','Normalized') 
xlabelC Actuator Thickness (% of Beam Thickness)') 
else 

plot(tp* 1 000, V_tlimit *tp. *K_Vl,'b~',tp* 1 000, V_tlimit*tp. *K_V 1 2,'m-.', . . . 

tp* 1000,V_tlimit*tp.*K_V2,'r-',tp* 1000,V_tlimit*tp.*K_V22,'g.') 
text(0.15,0.6,[num2str(wp/w*100) '% Width Actuator/' num2str(tb(l)*1000) 

' mm Bond Thickness'],Units','Normalized') 
xlabel('Actuator Thickness (mm)') 
end 

ylabel('Curvature (rad/m)') 

title(['Beam Curvature at ' num2str(Volt_mil) ' Volts/mil']) 
legend('Single','l Stacked',' Sym. Pair','Sym. Stacked',- 1) 
text(0. 1 5,0.9,Piezonam,Units','Normalized') 

text(0.15,0.8,['Ep = ' num2str(Ep) ' N/m^2 d31 = ' num2str(d31) ' m/V], ... 
'Units', 'Normalized') 

text(0.15,0.7,[num2str(t(l)*1000) ' mm Thick ' Beammatl ' Beam'],'Units', ... 

'Normalized') 

grid 

figure 

if Plottype = 'w' 

plot(wp./w* 1 00,StrainV 1 * le6,'b— ', wp./w* 1 00,StrainV12* 1 e6,'m-.', . . . 

wp./w* 100, Strain V2* le6,'r-',wp./w* 100,StrainV22* le6,'g.') 
text(0.15,0.6,[num2str(tp*1000) ' mm Thick Actuator/' num2str(tb(l)*1000) . 

' mm Bond Thickness'],'Units','Normalized') 
xlabelC Actuator Width (%)') 
elseif Axistype = 'p' 

plot(tp./t* 1 00,StrainV 1*1 e6,'b— ',tp./t* 1 00,StrainV 12*1 e6,'m-.', . . . 

tp ./t * 1 00, Strain V2* 1 e6,'r-',tp ./t * 1 00, Strain V22 * 1 e6,'g. ') 
text(0. 15,0. 6,[num2str(wp/w* 100) '% Width Actuator/' num2str(tb(l)*1000) 

' mm Bond Thickness'],'Units', 'Normalized') 
xlabelC Actuator Thickness (% of Beam Thickness)') 
else 

plot(tp * 1 000, Strain V 1*1 e6,'b— ',tp* 1 000, StrainV 1 2* 1 e6,'m-. ', . . . 

tp* 1000,StrainV2* le6,'r-',tp* 1000,StrainV22* le6,'g.') 
text(0. 15,0.6,[num2str(wp/w* 100) '% Width Actuator/" num2str(tb(l)*1000) 

' mm Bond Thickness'],'Units','Normalized') 
xlabel('Actuator Thickness (mm)') 
end 

ylabel('Strain/Voltage (microstrain/V)') 
title('Surface Axial Strain per Unit Actuator Voltage') 
legend('Single','l Stacked','Sym. Pair',' Sym. Stacked',- 1) 
text(0. 1 5,0.9,Piezonam,'Units','Normalized') 

text(0.15,0.8,['Ep = ' num2str(Ep) ' N/m'^2 d31 = ' num2str(d31) ' m/V], ... 
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' UnitsVNormalized') 

text(0.15,0.7,[num2str(t(l)*1000) ' mm Thick ' Beammatl ' Beam'],'Units', ... 

'Normalized') 

grid 

figure 

ifPlottype = 'w' 

plot(wp./w* 1 00,y0_l * 1 000,'b--',wp./w* 1 00,y0_l 2* 1 000,'m-.') 
text(0.15,0.6,[num2str(tp*1000) ' mm Thick Actuator/' num2str(tb(l)*1000) ... 

' mm Bond Thickness'],'Units','Normalized') 
xlabelC Actuator Width (%)') 
ylabel('yO (mm)') 
elseif Axistype = 'p' 

pIot(tp./t* 100,y0_l ./t* 100,'b-',tp./t* 1 00,y0_12./t* 100,'m-.') 
text(0.15,0.6,[num2str(wp/w*100) '% Width Actuator/' num2str(tb(l)*1000) ... 

' mm Bond Thickness'],Units','Normalized') 
xlabelC Actuator Thickness (% of Beam Thickness)') 
ylabel('yO (% of Beam Thickness)') 
else 

plot(tp* 1 000,y0_l * 1 000,'b-',tp* 1 000,y0_l 2* 1 000,'m-.') 
text(0.15,0.6,[num2str(wp/w*100) '% Width Actuator/' num2str(tb(l)*1000) ... 

' mm Bond Thickness'],'Units','Normalized') 
xlabelC Actuator Thickness (mm)') 
ylabel('yO (mm)') 
end 

titleCNeutral Axis Distance from Center of Beam') 

legend('Single','l Stacked',- 1) 

text(0. 1 5,0.9,Piezonam,Units','Normalized') 

text(0.15,0.8,['Ep = ' num2str(Ep) ' N/m'^2 d31 = ' num2str(d31) ' m/V'], ... 
'Units','Normalized') 

text(0.15,0.7,[num2str(t(l)*1000) ' mm Thick ' Beammatl ' Beam'],'Units', ... 

'Normalized') 

grid 
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% PROGRAM TO DETERMINE OPTIMUM ACTUATOR LOCATIONS 
% 

% This program determines the optimum actuator locations and input 
% voltages to achieve the minimum error between the actuator deformed 
% beam shape and a desired beam shape liinction for a cantilever 
% beam with piezoceramic actuators of specified length, width, thickness 
% and material properties bonded to both sides of the beam. 

% Set desired shape function and coefficients 
Shapfunc =1 % y = C*x^2 
% Shapfiinc = 2 % y = C*(l-cos(mx)) 

C = 0.001; 

% C = 0.0005; 
m = pi; 



% Set beam properties 
L = 1 ; % Beam length (m) 

w = 0.0508; % Beam width (m) 

% t = 0.001575; % Beam thickness (m) 

t = 0.001 575/2 % Beam thickness (m) 

E = 7.2el0; % Beam modulus of elasticity (N/m^2) 



% Set piezoceramic actuator properties 
% Navy Type II PZT 
n = 4 % Number of actuators 

la = 0. 1 *ones(n, 1); % Actuator lengths (m) 



wp = 0.0381; 
tp = 0.00027; 
Ep = 6.1el0; 
d31 = 1.71e-10, 
Vmax = 1 50; 
Vmin = -150; 



% Piezoceramic actuator width (m) 

% Piezoceramic actuator thickness (m) 

% Piezoceramic actuator modulus of elasticity (N/m^2) 
% Piezoelectric constant (m/V) 

% Maximum actuator input voltage (V) 

% Minimum actuator input voltage (V) 



% Enter properties of bond between piezoceramic actuator and beam 
wb = wp, % Bond width (m) (equal to actuator width) 

tb = 0; % Bond thickness (m) 

Eb = 1.78e9; % Bond modulus of elasticity (N/m^2) 



% Calculate properties of beam cross-section 

area = w.*t; % Beam cross-sectional area (m^2) 

areab = wb*tb; % Bond cross-sectional area (m^2) 

areap = wp*tp; % Piezoceramic actuator cross-sectional area (m^2) 

moi = w. *t.^3/ 1 2; % Beam cross-section moment of inertia (mM) 

% Bond cross-section moment of inertia (mM) 
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moib = wb*tb*(t^2/2 + t*tb + 2*tb^2/3), 

% Piezoceramic actuator cross-section moment of inertia (mM) 
moip = wp*tp*((t+2*tb)^2/2 + (t+2*tb)*tp + 2*tp''2/3); 

El = E*moi+Eb*moib+Ep*moip, % Beam stiffness 

K = Ep*d3 1 *wp*(t-t-2*tb+tp)/EI; % (Moment per unit voltage)/(Beam stiffness) 
K = K*ones(n,l); 

% Set options for unconstrained voltage optimization algorithm 
Options = foptions'; 

Options(14) = 2000; % Maximum number of iterations (0 -> 100*n) 

% Set initial actuator locations for optimization algorithm 
xO = [0.05 0.25 0.5 0.75]’ % 4 Actuators 

% Determine optimum voltage for initial actuator positions 
if Shapfunc = 1 

VO = fminsCbeamerVu',zeros(n,l),Options,[],xO,la,K,C,L,n) 

% VO = zeros(n,l) 
else 

VO = fmins('bercosVu',zeros(n,l),Options,[],xO,la,K,C,m,L,n) 
end 

% Calculate error cost function for initial actuator positions 
if Shapfimc = 1 

JO = beamerVu(V0,x0,la,K,C,L,n) 
else 

JO = bercosVu(V0,x0,la,K,C,m,L,n) 
end 

% Determine optimum actuator locations 
if Shapfunc == 1 

flops(O) % Reset floating point operations counter 

Optx = fmins('beamerxu',x0,Options,[],V0,la,K,C,L,n, Options) 

Ops = flops % Output number of floating point operations 
else 

flops(O) % Reset floating point operations counter 
% Unconstrained locations 

Optx = fmins('bercosxu',x0,Options,[],V0, la, K,C,m,L,n, Options) 

Ops = flops % Output number of floating point operations 
end 

% Determine optimum input voltages for optimum actuator locations 
if Shapfunc = 1 

OptV = fmins('beamerVu',VO,Options,[],Optx,la,K,C,L,n) 
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else 

OptV = fmins('bercosVu',VO,Options,[],Optx,la,K,C,m,L,n) 
end 

% Calculate error cost function for optimum actuator locations 
if Shapfunc = 1 

J = beamerVu(OptV,Optx,la,K,C,L,n) 
else 

J = bercosVu(OptV,Optx,la,K,C,m,L,n) 
end 

% Plot results 

% Set plot X range along beam length 
X = linspace(0,L,100); 

% Calculate beam displacement for optimum actuator locations and voltages 
y = beamdisp(x,Optx,OptV,la,n,K); 

% Calculate beam displacement at actuator actuator endpoints 
yact = beamdisp([Optx;Optx+la],Optx,OptV,la,n,K); 

% Calculate desired beam displacement along beam length 
if Shapfimc = 1 
yf= C*x.^2; 
else 

yf = C*(ones(size(x))-cos(m*x)); 
end 
figure 

plot(x,y,'r-',x,yf,'g.',[Optx;Optx+la],yact;r*') 

legend(' Actual VDesired',- 1 ) 

grid 

titleC Actual vs. Desired Beam Shapes') 
xlabeK'X (m)’) 

ylabel('Beam Displacement Y (m)') 
if n = 1 

text(0.1,0.9,['Optimum Values for 1 Actuator'],Units','Normalized') 
else 

text(0.1,0.9,['Optimum Values for ' int2str(n) ' Actuators'], 'Units', ... 

•Normalized') 

end 

if n = 1 

text(0.1,0.8,'['xl = ' num2str(0ptx(l)) ' m / VI = ' num2str(OptV(l)) ... 

' V],'Units','Normalized') 

text(0.1,0.7,['Error = ' num2str(J)],'Units', 'Normalized') 
elseif n = 2 

text(0.1,0.8,['xl = ' num2str(Optx(l)) ' m / VI = ' num2str(OptV(l)) ... 

' V],'Units','Normalized') 
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text(0. 1,0.7,['x2 = ' num2str(Optx(2)) ' m t Wl = ' num2str(OptV(2)) ... 

' V'],'Units','Nonnalized') 

text(0. l,0.6,['Error = ' num2str(J)],'Units','Normalized') 
elseif n = 3 

text(0.1,0.8,['xl = ' num2str(0ptx(l)) ' m / VI = ' num2str(0ptV(l)) ... 

' V'],Units','Normalized') 

text(0. 1,0.7,['x2 = ' num2str(Optx(2)) ' m / V2 = ' num2str(OptV(2)) ... 

' V'],'Units','Normalized') 

text(0.1,0.6,['x3 = ' num2str(Optx(3)) ' m / V3 = ' num2str(OptV(3)) ... 

' V'],'Units','Normalized') 

text(0. l,0.5,['Error = ' num2str(J)],TJnits','Normalized') 
elseif n = 4 

text(0.1,0.8,['xl = ■ num2str(0ptx(l)) ' m / VI = ' num2str(0ptV(l)) ... 

' V'],'Units','Normalized') 

text(0. 1,0.7,['x2 = ' num2str(Optx(2)) ' m / V2 = ' num2str(OptV(2)) ... 
' V],'Units','Normalized') 

text(0.1,0.6,['x3 = ' num2str(Optx(3)) ' m / V3 = ' num2str(OptV(3)) ... 
' V'],'Units','Normalized') 

text(0. 1,0. 5, ['x4 = ' num2str(Optx(4)) ' m / V4 = ' num2str(OptV(4)) ... 
' V],'Units','Normalized') 

text(0.1,0.4,['Error = ' num2str(J)],'Units',TSlormalized') 
else % (n = 5) 

text(0.1,0.8,['xl = ' nuin2str(Optx(l)) ' m / VI = ' num2str(OptV(l)) ... 
' V],'Units','Normalized') 

text(0.1,0.7,['x2 = ' num2str(Optx(2)) ' m / V2 = ' num2str(OptV(2)) ... 
' V'],'Units','Normalized') 

text(0.1,0.6,['x3 = ' num2str(Optx(3)) ' m / V3 = ' num2str(OptV(3)) ... 
' V'],Units','Normalized') 

text(0. 1,0.5, ['x4 = ' num2str(Optx(4)) ' m / V4 = ' num2str(OptV(4)) ... 
' V'],'Units','Normalized') 

text(0.1,0.4,['x5 = ' num2str(Optx(5)) ' m / V5 = ' num2str(OptV(5)) ... 
' V],’Units','Normalized') 

text(0.1,0.3,['Error = ' num2str(J)],'Units','Normalized') 
end 

% Allow data to be saved 
Savedata = input('Save data? (y/n) 
if Savedata = 'y' 

Savefile = input('Enter filename (without .MAT extension) 
eval(['save ' Savefile]) 
end 
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% FUNCTION TO DETERMINE ERROR BETWEEN ACTUAL BEAM SHAPE AND 
% DESIRED PARABOLIC SHAPE FOR UNCONSTRAINED ACTUATOR 
% LOCATION OPTIMIZATION 
% 

% function f = beamerxu(x,VO, la, K,C,L,n, Options) 

% 

% This function calculates the value of the cost function 'f 
% for a cantilever beam of length 'L' with 'n' attached 
% piezoceramic actuators at locations 'x' with lengths 'la' and 
% curvature per unit input voltage 'K' and a desired parabolic 
% beam shape function y = C*x''^2. 

function f = beamerxu(x,VO,la,K,C,L,n, Options) 

% Determine optimum actuator input voltages for input actuator 
% locations 'x' 

V = fmins('beamerVu',VO,Options,[],x,la,K,C,L,n); 

% Calculate summations for use in calculating beam shape 

% polynomial coeflRcients 

sumKVl = zeros(n,l); 

sumKVlxl = zeros(n,l); 

for i = l;n 

sumKVl(i) = sum(K(l;i).*V(l;i).*la(l:i)); % Sum ofKi*Vi*h 

sumKVlxl(i) = sum(K(l;i).*V(l;i).*la(l;i).*(x(l:i)+la(l:i)/2)); 

% Sum of Ki*Vi*li*(xi+li/2) 
end 

% Calculate beam shape polynomial coeflScients 
a2 = K.*V/2; 
al = -K.*V.*x; 

al(2:n) = al(2;n) + sumKVl(l:n-l); 

aO = K.*V.*x.^2/2; 

a0(2;n) = a0(2;n) - sumKVlxl(l:n-l); 

bl = sumKVl; 

bO = -sumKVlxl; 

% Create nxl vector containing value of 'C 
Cvect = C*ones(n,l); 

% Calculate value of cost function 
f = C^2*x(l)^5/5; 

f = f + sum((a2-Cvect).^2.*((x+la).'^5 - x.^5))/5; 
f = f + sum((a2-Cvect).*al.*((x+la).^4 - x.M))/2; 
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f= f + sum((2*(a2-Cvect).*aO+al.^2).*((x+la).''3 - x/'3))/3; 
f = f + sum(al .*a0.*((x+la).^2 - x.^2)); 
f = f + sum(a0.'^2.*la); 

f= f + C''2*sum(x(2;n).^5 - (x(l:n-l)+la(l;n-l)).^5)/5; 

f= f- C*sum(bl(l:n-l).*(x(2:n).^4 - (x(l;n-l)+la(l;n-l)).M))/2; 

f=f+ sum((-2*C*bO(l:n-l)+bl(l;n-l).''2).*(x(2;n).^3 - ... 

(x( 1 ; n- 1 )+la( 1 : n- 1 )) . ^^3 ))/3 ; 

f = f+ sum(bl(l:n-l).*b0(l;n-l).*(x(2:n).^2 - (x(l;n-l)+la(l.n-l)).^2)); 
f= f+ sum(b0(l;n-l).^2.*(x(2;n) - x(l:n-l)-la(l:n-l))); 
f = f + C^2*(L^5-(x(n)+la(n))^5)/5; 
f = f - C*bl(n)*(L^4-(x(n)+la(n))^4)/2; 
f = f + (-2*C*bO(n)+bl(n)^2)*(L^3-(x(n)+la(n))^3)/3; 
f = f + bl(n)*bO(n)*(L^2-(x(n)+la(n))^2); 
f = f + bO(n)^2*(L-x(n)-la(n)); 



% FUNCTION TO DETERMINE ERROR BETWEEN ACTUAL BEAM SHAPE AND 
% DESIRED PARABOLIC SHAPE FOR UNCONSTRAINED ACTUATOR 
% VOLTAGE OPTIMIZATION 
% 

% function f = BeamerVu(V,x,la,K,C,L,n) 

% 

% This function calculates the values of the cost function 'f 
% for a cantilever beam of 

% length U' with 'n' attached piezoceramic actuators and a 
% desired parabolic beam shape function y = C*x^2, where the 
% elements of the input vector 'V are the actuator input 
% voltages, the elements of 'x' are the distances 
% of the actuators from the root of the beam, and the 
% elements of 'la' are the actuator lengths. 

function f = BeamerVu(V,x,la,K,C,L,n) 

% x = Actuator locations (m) 

% V = Actuator voltages (V) 

% la = Actuator lengths (m) 

% Calculate summations for use in calculating beam shape 

% polynomial coefficients 

sumKVl = zeros(n,l); 

sumKVlxl = zeros(n,l); 

for i = 1 :n 

sumKVl(i) = sum(K(l;i).*V(l;i).*la(l:i)); % Sum of Ki*Vi*li 

sumKVlxl(i) = sum(K(l;i).*V(l:i).*la(l:i).*(x(l;i)+la(l:i)/2)); 
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end 



% Sum of Ki*Vi*li*(xi+li/2) 



% Calculate beam shape polynomial coefficients 
a2 = K.*V/2; 
al =-K.*V.*x; 

al(2;n) = al(2:n) + sumKVl(l;n-l); 

aO = K.*V.*x.^2/2; 

a0(2:n) = a0(2;n) - sumKVlxl(l;n-l); 

bl = sumKVl; 

bO = -sumKVbd; 

% Create nxl vector containing value of 'C 
Cvect = C*ones(n,l); 

% Calculate value of cost function 
f = C^2*x(l)^5/5; 

f = f + sum((a2-Cvect).^2.*((x+la).^5 - x.^5))/5; 
f = f + sum((a2-Cvect).*al.*((x+la).^4 - x.M))/2; 
f = f + sum((2*(a2-Cvect).*aO+al.^2).*((x+la).''3 - x.'^3))/3; 
f = f + sum(al.*a0.*((x+la).''2 - x.^2)); 
f = f + sum(a0.^2.*la); 

f = f+ C^2*sum(x(2:n).^5 - (x(l:n-l)+la(l:n-l)).^5)/5; 
f = f- C*sum(bl(l;n-l).*(x(2;n).M - (x(l;n-l)+la(l:n-l)).M))/2, 
f = f + sum((-2*C*bO(l;n-l)+bl(l;n-l).^2).*(x(2;n).^3 - ... 
(x(l:n-l)+la(l:n-l)).^3))/3; 

f=f+ sum(bl(l;n-l).*b0(l;n-l).*(x(2:n).^2 - (x(l:n-l)+la(l:n-l)).^2)); 
f= f + sum(b0(l;n-l).^2.*(x(2:n) - x(l:n-l)-la(l;n-l))); 
f = f + C^2*(L^5-(x(n)+la(n))^5)/5; 
f = f - C*bl(n)*(L^4-(x(n)+la(n)r4)/2; 
f = f + (-2*C*bO(n)+bl(n)^2)*(L^3-(x(n)+la(n))^3)/3; 
f = f + bl(n)*bO(n)*(L^2-(x(n)+la(n)r2); 
f = f + bO(n)^2*(L-x(n)-la(n)); 
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% FUNCTION TO DETERMINE ERROR BETWEEN ACTUAL BEAM SHAPE AND 
% DESIRED C*(I-cos(mx)) SHAPE FOR UNCONSTRAINED ACTUATOR 
% LOCATION OPTIMIZATION 
% 

% function f = bercosxu(x, VO,la,K,C,m,L,n,Options) 

% 

% This function calculates the value of the cost function 'f 
% for a cantilever beam of length 'L' with 'n' attached 
% piezoceramic actuators at locations 'x' with lengths 'la' and 
% curvature per unit input voltage 'K' and a desired sinusoidal 
% beam shape function y = C*(l-cos(mx)). 

function f = bercosxu(x,VO,la,K,C,m,L,n, Options) 

% Determine optimum actuator input voltages for input actuator 
% locations 'x' 

V = fitunsCbercosVu',VO,Options,[],x,la,K,C,m,L,n); 

% Calculate summations for use in calculating beam shape 

% polynomial coefficients 

sumKVl = zeros(n,l); 

sumKVTxl - zeros(n,l); 

for i = l:n 

sumKVl(i) = sum(K(l;i).*V(l:i).*la(l;i)); % Sum of Ki*Vi*li 

sumKVM(i) = sum(K(l.i).*V(l:i).*la(l:i).*(x(l:i)+la(l;i)/2)); 

% Sum of Ki*Vi*li*(xi+li/2) 
end 

% Calculate beam shape polynomial coefficients 
a2 = K.*V/2, 
al =-K.*V.*x; 

al(2:n) = al(2;n) + sumKVl(l;n-l); 

aO = K.*V.*x.^2/2; 

a0(2;n) = a0(2;n) - sumKVlxl(l;n-l); 

bl = sumKVl; 

bO = -sumKVlxl; 

% Create nxl vector containing value of 'C 
Cvect = C*ones(n,l); 

% Calculate value of cost function 

f = C^2*(x(l)-2/m*sin(m*x(l)) + L/2 + l/4/m*sin(2*m*L)); 
f = f + sum(ji2.^2.*((x+la).^5 - x.^5))/5; 
f = f+ sum(a2.*al.*((x+la).^4 - x.^4))/2; 
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f = f + sum((2*a2.*(aO-Cvect)+al.^2).*((x+la).^3 - x.^3))/3; 
f = f + sum(al.*(aO-Cvect).*((x+la).'^2 - x.^2)), 
f = f + sum((aO-Cvect).'^2.*la); 

f = f + 2*C*sum(a2.*(2/m^2*((x+la).*cos(m*(x+la))-x.*cos(m*x)))), 
f = f + 2*C*sum(a2.*(((x+la).^2/m-2/m^3*ones(n,l)).*sin(m*(x+la))- ... 

(x.''2/m-2/m^3 *ones(n, 1 )).*sin(m*x))); 
f = f + 2*C*sum(al.*((cos(m*(x+la))-cos(m*x))/m^2 + ... 

((x+la). *sin(m*(x-Ha))-x. *sin(m*x))/m)); 
f = f + 2/m*C*sum((aO-Cvect).*(sin(m*(x+la))-sin(m*x))); 
f=f+ sum(bl(l;n-l).'^2.*(x(2;n).^3 - (x(l:n-l)+la(l:n-l)).^3))/3; 
f=f + sum(bl(l:n-l).*(bO(l:n-l)-Cvect(l:n-l)).*(x(2:n).'^2 - ... 

(x( 1 : n- 1 )+la( 1 : n- 1 )). ^^2)); 

f = f + sum((bO(l:n-l)-Cvect(l:n-l)).^2.*(x(2;n) - x(l:n-l)-la(l:n-l))); 
f = f + 2*C*sum(bl(l;n-l).*(cos(m*x(2:n))-cos(m*(x(l:n-l)+la(l:n-l)))))/m^2; 
f = f + 2*C*sum(bl(l;n-l).*(x(2:n).*sin(m*x(2;n))-(x(l;n-l)+la(l:n-l)).* ... 

sin(m*(x(l ;n-l )+la(l :n-l)))))/m; 
f= f + 2*C/m*sum((bO(l:n-l)-Cvect(l;n-l)).*(sin(m*x(2:n))- ... 

sin(m* (x( 1 ; n- 1 )+la( 1 ; n- 1 ))))); 
f = f + bl(n).^2.*(L.^3 - (x(n)+la(n)).^3)/3; 
f = f + bl(n).*(bO(n)-C).*(L.^2 - (x(n)+la(n)).^2); 
f = f + (bO(n)-C).'^2.*(L - x(n)-la(n)); 
f = f + 2*C*bl(n).*(cos(m*L)-cos(m*(x(n)+la(n))))/m^2; 
f = f + 2*C*bl(n).*(L.*sin(m*L)-(x(n)+la(n)).*sin(ni*(x(n)+la(n))))/m; 
f = f + 2*C/m*(bO(n)-C).*(sin(m*L)-sin(m*(x(n)+la(n)))); 



% FUNCTION TO DETERMINE ERROR BETWEEN ACTUAL BEAM SHAPE AND 
% DESIRED C*(l-cos(mx)) SHAPE FOR UNCONSTRAINED ACTUATOR 
% VOLTAGE OPTIMIZATION 
% 

% function f = bercosVu(x,VO,la,K,C,m,L,n, Options) 

% 

% This function calculates the value of the cost function 'f 
% for a cantilever beam of length 'L' with 'n' attached 
% piezoceramic actuators at locations 'x' with lengths 'la' and 
% curvature per unit input voltage 'K' and a desired sinusoidal 
% beam shape function y = C*(l-cos(mx)). 

function f = bercosVu(V,x,la,K,C,m,L,n,Options) 



% X = Actuator locations (m) 
% V = Actuator voltages (V) 
% la = Actuator lengths (m) 
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% Calculate summations for use in calculating beam shape 

% polynomial coefficients 

sumKVl = zeros(n,l); 

sumKVlxl = zeros(n,l); 

for i = l:n 

sumKVl(i) = sum(K(l:i).*V(l:i).*la(l:i)); % Sum ofKi*Vi*li 

sumKVlxl(i) = sum(K(l;i).*V(l:i).*la(l;i).*(x(l;i)+la(l;i)/2)); 

% Sum of Ki*Vi*li*(xi+li/2) 
end 

% Calculate beam shape polynomial coefficients 
a2 = K.*V/2; 
al =-K.*V.*x; 

al(2;n) = al(2:n) + sumKVl(l:n-l); 

aO = K.*V.*x.'^2/2; 

a0(2:n) = a0(2;n) - sumKVlxl(l;n-l); 

bl = sumKVl; 

bO = -sumKVlxl; 

% Create nxl vector containing value of 'C 
Cvect = C*ones(n,l); 

% Calculate value of cost function 

f = C^2*(x(l)-2/m*sin(m*x(l)) + L/2 + l/4/m*sin(2*m*L)); 
f = f + sum(a2.^2.*((x+la).^5 - x.^5))/5; 
f = f + sum(a2.*al.*((x+la).M - x.M))/2; 
f = f + sum((2*a2.*(aO-Cvect)+al.'^2).*((x+la).''3 - x.'^3))/3; 
f = f + sum(al .*(aO-Cvect).*((x+la).^2 - x.^2)); 
f = f + sum((aO-Cvect).'^2.*la); 

f = f + 2*C*sum(a2.*(2/m^2*((x+la).*cos(m*(x+la))-x.*cos(m*x)))); 
f = f + 2*C*sum(a2.*(((x+la).^2/m-2/m'^3*ones(n,l)).*sin(m*(x+la))- ... 

(x.^2/m-2/m^3 * ones(n, 1 )). *sin(m*x))); 
f = f + 2*C*sum(al.*((cos(m*(x+la))-cos(m*x))/m^2 + ... 

((x+la). * sin(m*(x+la))-x. *sin(m*x))/m)); 
f = f + 2/m*C*sum((aO-Cvect).*(sin(m*(x+la))-sin(m*x))); 
f= f+ sum(bl(l;n-l).^2.*(x(2:n).^3 - (x(l;n-l)+la(l:n-l)).^3))/3; 
f= f+ sum(bl(l:n-l).*(bO(l:n-l)-Cvect(l;n-l)).*(x(2:n).^2 - ... 
(x(l;n-l)+la(l:n-l)).^2)); 

f= f+ sum((bO(l:n-l)-Cvect(l;n-l)).^2.*(x(2:n) - x(l;n-l)-la(l:n-l))); 

f= f + 2*C*sum(bl(l:n-l).*(cos(m*x(2:n))-cos(m*(x(l;n-l)+la(l:n-l)))))/m^2; 

f = f + 2*C*sum(bl(l:n-l).*(x(2:n).*sin(m*x(2;n))-(x(l;n-l)+la(l:n-l)).* ... 

sin(m* (x( 1 ; n- 1 )+la( 1 ; n- 1 )))))/m; 
f = f+ 2*C/m*sum((bO(l;n-l)-Cvect(l:n-l)).*(sin(m*x(2:n))- ... 
sin(m*(x( 1 : n- 1 )+la( 1 : n- 1 ))))); 
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f = f + bl(n).^2.*(L.'^3 - (x(n)+la(n)).^3)/3; 

f = f + bl(n).*(bO(n)-C).*(L.^2 - (x(n)+la(n)).^2), 

f = f + (bO(n)-C).^2.*(L - x(n)-la(n)); 

f = f + 2*C*bl(n).*(cos(m*L)-cos(m*(x(n)+la(n))))/m'^2; 

f = f + 2*C*bl(n).*(L.*sin(m*L)-(x(n)+la(n)).*sin(m*(x(n)+la(n))))/m; 

f = f + 2*C/m*(bO(n)-C).*(sin(m*L)-sin(m*(x(n)+la(n)))); 



% FUNCTION TO DETERMINE THE DEFLECTED SHAPE OF A CANTILEVER 
% BEAM DUE TO ATTACHED PIEZOCERAMIC ACTUATORS 
% 

% function y = beamdisp(x,xa,v,la,n,K) 

% 

% This function returns a vector of transverse displacements, y, 

% of a cantilever beam with 'n' attached piezoceramic actuators, 

% given the starting locations of the actuators, 'xa', the lengths 
% of the actuators, 'la', the actuator input voltages, 'v', and the 
% vector 'K' containing the actuator curvatures per unit voltage. 

% 'x' is a vector containing a range of axial distances from the 
% root of the beam at which 'y' is to be evaluated and 'xa', 'v' 

% and 'la' are all nxl vectors. 

function y = beamdisp(x,xa,v,la,n,K) 
y = zeros(size(x)); 
for m = 1 :length(x) 

if x(m) <= xa(l ) % x prior to start of actuator 

y(m) = 0; 

elseif xa(l) < x(m) & x(m) <= xa(l)+la(l) 
y(m) = K(l)*v(l)/2*(x(m)-xa(l))'^2; 
elseif x(m) <= xa(n)+la(n) 
for k= l:n 

if xa(k) <= x(m) & x(m) <= xa(k)+la(k) 
y(m) = K(k)*v(k)/2*(x(m)-xa(k))''2+sum(K(l;k-l).*v(l:k-l).* ... 
la( 1 ; k- 1 ). *(x(m) *ones(size(xa( 1 ; k- 1 )))-xa( 1 ; k- 1 )-la( 1 ; k- 1 )/2)), 
elseif k<n 

if xa(k)+la(k) <= x(m) & x(m) <= xa(k+l) 
y(m) = sum(K(l:k).*v(l;k).*la(l;k).*(x(m)*ones(size(xa(l:k)))- ... 
xa(l:k)-la(l:k)/2)), 
end 
end 
end 
else 

y(m) = sum(K(l:n).*v(l;n).*la(l;n).*(x(m)*ones(size(xa(l:n)))- ... 
xa( 1 : n)-la( 1 ; n)/2)); 
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end 

end 



% PROGRAM TO DETERMINE OPTIMUM ACTUATOR LOCATIONS 
% FROM A GIVEN DISCRETE SET OF POSSIBLE LOCATIONS 
% 

% This program determines the optimum actuator location and input 
% voltage from a specified discrete range of possible actuator 
% locations to achieve the minimum error between the actuator deformed 
% beam shape and a desired parabolic shape function for a cantilever 
% beam with piezoceramic actuators bonded to both sides of the beam. 

% 1, 2 or 3 actuators can be specified 

% Coefficient of desired parabolic shape function 
C = 0.001; 



% Enter beam properties 
L = 1 ; % Beam length (m) 

w = 0.0508; % Beam width (m) 

t = 0.001575; % Beam thickness (m) 

E = 7.2el0; % Beam modulus of elasticity (N/m^2) 



% Enter piezoceramic actuator properties 
% Navy Type II PZT 



n = 3; % Number of actuators 

la = 0 . 1 *ones(n, 1 ); % Actuator lengths (m) 

% Piezoceramic actuator width (m) 

% Piezoceramic actuator thickness (m) 

% Piezoceramic actuator modulus of elasticity (N/m^2) 
% Piezoelectric constant (m/V) 



wp = 0.0381; 
tp = 0.00027; 
Ep = 6.1el0; 
d31 = 1.71e-10; 



% Enter properties of bond between piezoceramic actuator and beam 
wb = wp; % Bond width (m) (equal to actuator width) 

tb = 0; % Bond thickness (m) 

Eb = 1 .78e9; % Bond modulus of elasticity (N/m^2) 



% Calculate properties of beam cross-section 

area = w.*t; % Beam cross-sectional area (m^2) 

areab = wb*tb; % Bond cross-sectional area (m^2) 

areap = wp*tp; % Piezoceramic actuator cross-sectional area (m^2) 

moi = w. *t.^3/l 2; % Beam cross-section moment of inertia (m^4) 

% Bond cross-section moment of inertia (m'^4) 
moib = wb*tb*(t''2/2 + t*tb + 2*tb"'2/3); 
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% Piezoceramic actuator cross-section moment of inertia (mM) 
moip = wp*tp*((t+2*tb)^2/2 + (t+2*tb)*tp + 2=^tp^2/3); 

El = E*moi+Eb*moib+Ep*moip; % Beam stifihess 

K = Ep*d3 1 *wp*(t+2*tb+tp)/EI; % (Moment per unit voltage)/(Beam stiffness) 
K = K*ones(n,l); 

% Set run number, minimum and maximum locations for actuators 
% and spacing between actuator locations in range of positions 
% evaluated between minimum and maximum locations for each actuator 
Run = 6; % Run number 

Savefile = 'NTn3_6'; 

xmin = [0.055 0.355 0.655]; % Run 6 (0.005 spacing) 

xmax = [0.095 0.395 0.695]; % Run 6 (0.005 spacing) 

Spacing = [0.005 0.005 0.005]; % Run 6 

% Set initial voltages for unconstrained voltage optimization algorithm 
VO = zeros(n,l); 

% Set options for unconstrained voltage optimization algorithm 
Options = foptions'; 

Options(14) = 2000; % Maximum number of iterations (0 -> 100*n) 

% Divide each actuator location range into up to 9 equally spaced 
% possible locations separated by the interval set in 'Spacing' 
if n = 1 % 1 Actuator 

xrange = [xmin:Spacing:xmax]; 
elseif n = 2 % 2 Actuators 

xrange = [xmin(l):Spacing(l):xmax(l);xmin(2):Spacing(2):xmax(2)], 
else % 3 Actuators 

xrange = [xmin(l);Spacing(l):xmax(l);xinin(2):Spacing(2):xmax(2); ... 
xmin(3); Spacing(3):xmax(3)]; 

end 

% Evaluate optimum actuator voltages and beam shape error 
% for each combination of actuator locations 
X = zeros(n,l); 

Maxindex = length(xrange); 
if n = 1 % 1 Actuator 

Ermatrix = zeros(Maxindex,l); 

Voltage 1 = zeros(Maxindex,l); 
for k = l;Maxindex 
fprintf(['\nk = ' int2str(k) '\n']) 

X = xrange(k); 

V = fmins('beamerVu',V0,Options,[],x,la,K,C,L,n); 
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Voltagel(k) = V; 

Ermatrix(k) = beamerVu(V,x,la,K,C,L,n); 
ifk--= 1 

if Ennatrix(k) < Minerror 
Optloc = k; 

OptV= V; 

Minerror = Ermatrix(k) 
end 
else 

Optloc = 1 ; 

OptV = V; 

Minerror = Ermatrix(k) 
end 
end 
clear k 

Optx - xrange(Optloc) 

OptV 

Minerror 

elseif n — 2 % 2 Actuators 

Emiatrix = zeros(lO); 

Voltagel = zeros(lO); 

Voltage2 = zeros(lO); 

Voltsum = zeros(lO); 
for j = IrMaxindex 
x(l) = xrange(lj); 
for k = liMaxindex 

fprintf(['\njk = ' int2str(j) ' ' int2str(k) '\n']) 
x(2) = xrange(2,k); 
if (x(l)+la(l)<=x(2)) 

V = finins('beamerVu',VO,Options,[],x,la,K,C,L,n); 
Voltagel(j,k) = V(l); 

Voltage2(j,k) = V(2); 

Ermatrix(j,k) = beamerVu(V,x,la,K,C,L,n); 
ifj+k ~= 2 

if Ermatrix(j,k) < Minerror 
Optloc = [j k]; 

OptV = V; 

Minerror = Ermatrix(j,k) 
end 
else 

Optloc =[1 1]; 

OptV = V; 

Minerror = Ermatrix(j,k) 
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end 

else 

Voltage l(j,k) = 0; 

Voltage2(j,k) = 0; 

Ermatrix(j,k) = 0; 
end 
end 
end 

clear j k 

Optx = [xrange(l,Optloc(l));xrange(2,Optloc(2))] 

OptV 

Minerror 

Voltsum = Voltagel+Voltage2; 

else % 3 Actuators 

Ermatrix = zeros(100,10); 

Voltagel =zeros(100,10); 

Voltage2 = zeros(100,10); 

Voltages =zeros(100,10); 

Voltsum = zeros(100,10); 
for i = EMaxindex 
x(l) = xrange(l,i); 
for j = EMaxindex 
x(2) = xrange(2,j); 
for k = 1 ;Maxindex 

^rintf(['\nijk = ' int2str(i) ' ' int2str(j) ' ' int2str(k) '\n']) 
x(3) - xrange(3,k); 

if (x(l)+la(l)<=x(2)) & (x(2)+la(2)<=x(3)) 

V = finins('beamerVu',VO,Options,[],x,la,K,C,L,n); 
Voltagel (10*(i-l)+j,k) = V(l); 
Voltage2(10*(i-l)+j,k) = V(2); 
Voltage3(10*(i-l)+j,k) = V(3); 
Ermatrix(10*(i-l)+j,k) = beamerVu(V,x,la,K,C,L,n); 
if i+j+k ~= 3 

if Ermatrix(10*(i-l)+j,k) < Minerror 
Optloc = [i j k]; 

OptV = V; 

Minerror = Ermatrix(10*(i-l)+j,k) 
end 
else 

Optloc =[11 1]; 

OptV = V; 

Minerror = Ermatrix(10*(i-l)+j,k) 
end 
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else 

Voltagel(10*(i-l)+j,k) = 0; 

Voltage2(10*(i-l)+j,k) = 0; 

Voltage3(I0*(i-l)+j,k) = 0; 

Ermatrix(10*(i-l)+j,k) = 0; 
end 
end 
end 
end 

clear i j k 

Optx = [xrange( 1 ,Optloc( 1 ));xrange(2,Optloc(2));xrange(3,Optloc(3))] 

OptV 

Minerror 

Voltsum == Voltage l+Voltage2+Voltage3; 
end 

% Plot results 
X = linspace(0,L,100); 

y = beamdisp(x, Optx, Opt V,la,n,K); % Fixed length actuators 

yact = beamdisp([Optx;Optx+la],Optx,OptV,la,n,K); 
yf = C*x.^2; 
figure 

plot(x,y,'r-’,x,yf,'g.',[Optx;Optx+la],yact,'r*') 

legend ('ActualVDesired',- 1 ) 

grid 

titleC Actual vs. Desired Beam Shapes') 
xlabel('X (m)') 

ylabel('Beam Displacement Y (m)') 
if n = 1 

text(0.1,0.9,['Run ' int2str(Run) ' Optimum Values for 1 Actuator;'], ... 
'Units','Normalized') 

text(0. 1,0.8, ['xl = ' num2str(Optx(l)) ' m / VI = ' num2str(OptV(l)) ... 
' V'],'Units', 'Normalized') 

text(0.1,0.7,['Error = ' num2str(Minerror,4)],'Units','Normalized') 
elseif n = 2 

text(0.1,0.9,['Run ' int2str(Run) ' Optimum Values for 2 Actuators;'], ... 
'Units','Normalized') 

text(0.1,0.8,['xl = ' num2str(Optx(l)) ' m / VI = ' num2str(OptV(l)) ... 
' V],'Units','Normalized') 

text(0.1,0.7,['x2 = ' num2str(Optx(2)) ' m / VI = ' num2str(OptV(2)) ... 
' V'],'Units','Normalized') 

text(0. 1,0.6, ['Error = ' num2str(Minerror,4)],Units','Normalized') 
elseif n = 3 

text(0.1,0.9,[Run ' int2str(Run) ' Optimum Values for 3 Actuators;'], ... 
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'Units','Nonnalized') 

text(0.1,0.8,['xl = ' num2str(0ptx(l)) ' m / VI = ' num2str(0ptV(l)) ... 
' V'],'Units','Normalized') 

text(0.1,0.7,['x2 = ' num2str(Optx(2)) ' m / VI = ' num2str(OptV(2)) ... 
' V'],'Units',TSfornialized') 

text(0.1,0.6,['x3 = ' num2str(Optx(3)) ' m / VI = ' num2str(OptV(3)) ... 
' V'],'Units','Normalized') 

text(0.1,0.5,['Error = ' num2str(Minerror,4)],'Units','Normalized') 
end 

% Save data to .MAT file 
eval(['save ' Savefile]) 



% PROGRAM TO DETERMINE OPTIMUM ACTUATOR INPUT VOLTAGES 
% FROM A GIVEN DISCRETE SET OF POSSIBLE VOLTAGES 
% 

% This program determines the optimum actuator input voltages 
% from a specified discrete range of possible actuator 
% voltages to achieve the minimum error between the actuator deformed 
% beam shape and a desired parabolic shape function for a cantilever 
% beam with piezoceramic actuators bonded to both sides of the beam. 

% 1, 2 or 3 actuators can be specified 

% Coefficient of desired parabolic shape function 
C = 0.001; 



% Enter beam properties 
L = 1 ; % Beam length (m) 

w = 0.0508; % Beam width (m) 

t = 0.001575; % Beam thickness (m) 

E = 7.2el0; % Beam modulus of elasticity (N/m^2) 



% Enter piezoceramic actuator properties 
% Navy Type II PZT 
n = 3; % Number of actuators 

la = 0. 1 *ones(n, 1 ); % Actuator lengths (m) 

% Piezoceramic actuator width (m) 

% Piezoceramic actuator thickness (m) 

% Piezoceramic actuator modulus of elasticity (N/m^2) 
% Piezoelectric constant (nW) 

% Maximum actuator input voltage (V) 

% Minimum actuator input voltage (V) 



wp = 0.0381; 
tp = 0.00027; 
Ep = 6.1el0; 
d31 = 1.71e-10; 
Vmax = 150; 
Vmin = -150, 
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% Enter properties of bond between piezoceramic actuator and beam 
wb = wp; % Bond width (m) (equal to actuator width) 

tb = 0, % Bond thickness (m) 

Eb = 1.78e9; % Bond modulus of elasticity (N/m^2) 



% Calculate properties of beam cross-section 

area = w.*t; % Beam cross-sectional area (m^2) 

areab = wb*tb; % Bond cross-sectional area (m^2) 

areap = wp*tp; % Piezoceramic actuator cross-sectional area (m^2) 

moi = w.*t.^3/12; % Beam cross-section moment of inertia (mM) 

% Bond cross-section moment of inertia (mM) 
moib = wb*tb*(t^2/2 + t*tb + 2*tb^2/3); 

% Piezoceramic actuator cross-section moment of inertia (mM) 
moip = wp*tp*((t+2*tb)''2/2 + (t+2*tb)*tp + 2*tp^2/3), 

El = E*moi+Eb*moib+Ep*moip; % Beam stiffhess 

K = Ep*d3 1 * wp*(t+2*tb+tp)/EI; % (Moment per unit voltage)/(Beam stiffiiess) 
K = K*ones(n,l); 



% Set run number, minimum and maximum voltages for actuators 
% and increment between actuator voltages in range of positions 
% evaluated between minimum and maximum locations for each actuator 
Run =1; % Run number 

Savefile = 'NTU1_1 V; 

Vmin = [10]; % Run 1 (5 V increment) 

Vmax = [50]; % Run 1 (5 V increment) 

Incremnt = [5]; % Run 1 

% Set actuator locations 
ifn = 1 
X = 0.2420; 
elseif n = 2 
x = [0.1220;0.5394]; 
else 

x = [0.0728;0.3678;0.6641]; 
end 



% Set initial voltages for unconstrained voltage optimization algorithm 
VO = zeros(n,l); 

% Set options for unconstrained voltage optimization algorithm 
Options = foptions'; 

Options(14) = 2000; % Maximum number of iterations (0 -> 100*n) 

% Divide each actuator voltage range into up to 9 equally spaced 
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% possible locations separated by the interval set in 'Spacing' 
if n = 1 % 1 Actuator 

Vrange = [Vminiincremnt; Vmax]; 
elseif n = 2 % 2 Actuators 

Vrange = [Vmin(l):Incremnt(l): Vmax(l);Vmin(2);Incremnt(2); Vmax(2)], 
else % 3 Actuators 

Vrange = 

[Vmin(l):Incremnt(l):Vmax(l);Vmin(2);Increntnt(2);Vmax(2);Vniin(3):Incremnt(3):Vma 

x(3)]; 

end 

% Evaluate optimum actuator voltages and beam shape error 
% for each combination of actuator locations 
Maxindex = length( Vrange); 
if n = 1 % 1 Actuator 

Ermatrix = zeros(Maxindex,l); 

Voltage 1 = zeros(Maxindex,l); 
for k = 1 : Maxindex 
fprintf(['\nk = ' int2str(k) '\n']) 

V = Vrange(k); 

Ermatrix(k) = beamerVu(V,x,la,K,C,L,n); 
if k “ 1 

if Ermatrix(k) < Minerror 
OptV = V; 

Minerror = Ermatrix(k) 
end 
else 

OptV = V; 

Minerror = Ermatrix(k) 
end 
end 
clear k 
OptV 
Minerror 

elseif n = 2 % 2 Actuators 

Ermatrix = zeros(lO); 
for j = 1 ;Maxindex 
for k = 1 :Maxindex 

fprintf(['\njk = ’ int2str(j) ' ' int2str(k) '\n']) 

V = [Vrange(l j);Vrange(2,k)]; 

Ermatrix(j,k) = beamerVu(V,x,la,K,C,L,n); 
if j+k ~ 2 

if Ermatrix(j,k) < Minerror 
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OptV = V; 

Minerror = Ermatrix(j,k) 
end 
else 

OptV = V; 

Minerror = Ermatrix(j,k) 
end 
end 
end 

clear j k 

OptV 

Minerror 

else % 3 Actuators 

Ermatrix = zeros(100,10), 
for i = 1 ;Maxindex 
for j = l:Maxindex 
for k = 1 ;Maxindex 

fprintf([’\nijk = ' int2str(i) ' ' int2str(j) ’ ' int2str(k) V]) 
V = [Vrange(l,i);Vrange(2,j);Vrange(3,k)]; 
Ermatrix(10*(i-l)+j,k) = beamerVu(V,x,la,K,C,L,n); 
ifi+j+k-= 3 

if Ermatrix(10*(i-l)+j,k) < Minerror 
OptV = V; 

Optindex = [i j k]; 

Minerror = Ermatrix(10*(i-l)+j,k) 
end 
else 

OptV = V; 

Optindex = [i j k]; 

Minerror = Ermatrix(10*(i-l)+j,k) 
end 
end 
end 
end 

clear i j k 
OptV 
Minerror 
end 

% Plot results 

if n = 1 % 1 Actuator 

% Plot variation in shape error vs. actuator voltage 
figure 
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plot(Vrange,Ermatrix) 

title('Error vs. Actuator Input Voltage for 1 Actuator') 

xlabelC Actuator Voltage (V)') 

ylabel('Error') 

text(0.15,0.95,['Run ' int2str(Run)],'Units','Normalized') 

text(0.1 5,0. 85,[' Actuator Location: ' num2str(x) ' m'], 'Units', 'Normalized') 

grid 

elseif n — 2 % 2 Actuators 

% Plot variation in shape error vs. actuator positions 
figure 

meshz(Vrange(2,:), Vrange( 1 ,;),Ermatrix(l iMaxindex, 1 .Maxindex)) 
title('Error vs. Actuator Input Voltages for 2 Actuators') 
xlabel('Voltage 2 (V)') 
ylabel('Voltage 1 (V)') 

text(0.15,0.95,['Run ' int2str(Run)],'Units', 'Normalized') 
text(0. 1 5,0 . 8 5, ['Actuator Locations: ' num2str(x( 1 )) ' m ' . . . 
num2str(x(2)) ' m'],'Units',Normalized') 

else % 3 Actuators 

% Plot variation in shape error for Actuator 1 fixed 

% at optimum location 

figure 

meshz( Vrange(3 , : ), Vrange(2, :),Ermatrix( 1 0 * (Optindex( I )- 1 )+ 1 : ... 

1 0 * (Optindex( 1 )- 1 )+Maxindex, 1 : Maxindex)) 
title('Error vs. Actuator 2 & 3 Voltage at Optimum Actuator 1 Voltage') 
xlabelC Voltage 3 (V)') 
ylabel('Voltage 2 (V)') 

text(0.15,0.95,['Run ' int2str(Run)],'Units','Normalized') 
text(0.15,0.85,[' Actuator Locations: ' num2str(x(l)) ' m ' ... 
num2str(x(2)) ' m ' num2str(x(3)) ' m'],'Units','Normalized') 

% Plot variation in shape error for Actuator 2 fixed 

% at optimum location 

figure 

meshz( Vrange(3 , : ), V range( 1 ,: ),Ermatrix( 1 0 * [0: Maxindex- 1 ]+ . . . 

Optindex(2)*ones( 1, Maxindex), 1 :Maxindex)) 
title('Error vs. Actuator 1 & 3 Voltage at Optimum Actuator 2 Voltage') 
xlabel('Voltage 3 (V)') 
ylabel('Voltage 1 (V)') 

text(0.15,0.95,['Run ' int2str(Run)],'Units','Normalized') 
text(0.15,0.85,['Actuator Locations: ' num2str(x(l)) ' m ' ... 
num2str(x(2)) ' m ' num2str(x(3)) ' m'],'Units','Normalized') 
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% Plot variation in shape error for Actuator 3 fixed 

% at optimum location 

figure 

Ermesh = reshape(Ermatrix(;,Optindex(3)),10,10)', 
meshz(Vrange(2,;),Vrange(l,;),Ermesh(l ;Maxindex,l ;Maxindex)) 
title('Error vs. Actuator 1 & 2 Voltage at Optimum Actuator 3 Voltage') 
xlabel('Voltage 2 (V)') 
ylabel('Voltage 1 (V)') 

text(0.15,0.95,['Run ' int2str(Run)],'Units','Normalized') 
text(0. 1 5, 0.85,[' Actuator Locations; ' num2str(x(l)) ' m ' 
num2str(x(2)) ' m ' num2str(x(3)) ' m'],'Units','Normalized') 
end 

% Save data to .MAT file 
eval(['save ' Savefile]) 
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